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1. Introduction. First we describe the channels we will discuss in this paper. The 
alphabet in which we send words contains just two letters, which we denote by 0, 1. 
Thus a sent word of length n is a sequence of n letters, each letter being a zero or a one. 
If we send the word (x, , --- , 2,), then the received word is (Y; , --- , Y,), 
where Y, . +++ , Y, are independent chance variables, the probability distribution of 
Y, depending only on the value of the parameter xz; . If there are just two possible 
values for Y; , which we can assume are 0, 1, the channel is called a “binary channel.”’ 
If Y, has a probability density function whether 2; is zero or one, the channel is called 
a ‘‘semi-continuous channel.”’ For a binary channel, f,,(y) denotes P(Y; = y, when the 
ith letter sent is x,). For a semi-continuous channel, f,,(y) denotes the probability density 
function of Y; when the 7th letter sent is z; . 

A binary channel is called “‘symmetric’”’ if f.(0) = f,(1) (and therefore f,(1) = f,(0)). 
A semi-continuous channel is called “symmetric” if the probability distribution of 
2[f.(1) + f,(¥)])'fo(Y) when Y has the probability density function f,(y) is the same 
as the probability distribution of 2[f.(Y) + f,(Y)]"'fi1(Y) when Y has the probability 
density function f,(y). 

A code of length L, word length n, and probability of error not greater than } is a 
sequence of L pairs (u, , Ay), «++ , (uz , Ax) with the following properties: 

(a) For each 7, u,; is a sequence of zeros and ones, n symbols altogether; 

b) For each 7, A, is a collection of words of length n, each being a possible received 
word in the particular problem under consideration; 

c) For each 7, when u; is the sent word, P((Y,, --- , Y,) isin A;) > 1— A”; 

d) The sets A, , --- , A, are disjoint. 

The use of such a code is well known. If the received word is in A; , the receiver 
assumes that the word u; was actually sent. Then, no matter which of the words u, , --- , 
u, is sent, the probability that the receiver will be in error is not greater than X. 

A ‘“‘converse to the coding theorem”’ is an inequality giving an upper bound for ZL, 
or, alternatively, for log.l. (In this paper, each log is to the base e, unless another base 
is explicitly indicated.) Since several types of converse have appeared in the literature, 
it seems worthwhile to classify them carefully. The following classification seems useful. 


‘Received May 25, 1959. This research was performed under contract with the Office of Naval 
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C denotes the capacity of the channel as defined in Feinstein [2]. A ‘‘weak converse” 
states the following. For any fixed « > 0, there is a positive value Z, such that there 
cannot exist a code of length 2"‘°*® and probability of error not greater than Z, for large 
n. A “strong converse” states the following. For any fixed « > 0, and any fixed X in the 
open interval (0, 1), there cannot exist a code of length 2"‘°*® and probability of error 
not greater than X, for n large. A “stronger converse” states the following. For any 
fixed \ in the open interval (0, 1), there is a constant K, such that there cannot exist a 
code of length 2°°*™"’” and probability of error not greater than X, for large n. 
Clearly, each of the converses listed is stronger than its predecessors. An interesting 
weak converse applicable to many types of channels is given in Theorem 5 of [3]. Wolfo- 
witz [4] has proved the stronger converse for a binary channel, with A, a positive constant, 
and [5] has proved the strong converse for a semi-continuous channel. It is the purpose 
of the present paper to show that for both binary symmetric channels and semi-continuous 
symmetric channels, the stronger converse can be proved with K, a negative constant 
if \ < 4. This makes the converse still stronger. 
We define 8, by the equation 


B 


[ (27)? exp (—30) dt = 1—X. 
Then if \ < 4, 8, > 0. 

2. The binary symmetric channel. 
f,(1) by q, and 1 — q by p. We assume that q is in the open interval (3, 1). The capacity 
C for this channel is 1 + p log.p + q logzg. 

We have the following theorem. For any 6 > 0, there is a number n, such that ifn > ns , 
then the length L of any code of word length n and probability of error not greater than X 
must satisfy the inequality log.L. < nC — n‘”* [8,(pq)'”” log.(q/p) — 4]. 


We denote the common value of f,(0) and 


Proof. If the word u; is sent, the most probable received word is u, itself, with prob- 


sie i n\ ia ae 
ability q”; the (") words differing from u; in exactly one letter are tied for next most 


; : — a n es ‘ 
probable received word, each having probability pq"~’; the (") words differing from 


u; in exactly two letters are tied for next most probable received word, each having 


n-2 


probability p’q”’; ete. 
We define K as the largest integer such that 


3 if] 
ct) pighi<s1—2. (1) 


=0 \J 


Denote by M; the number of words in A; . Then we must have 


K - 
M,> > ("), 
i=0 \J 
‘ 
since even the >> 


=0 


i . . 1 
( which are the most probably received words have a total prob- 
ability no greater than 1 — X. 

Since there are 2” different words of length n, and since A, , «++ , A, are disjoint, 


we have 
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I K = 
2> DM.>Ly ("), or L<2/ 2, (‘) ” r/'(x). 


i=0 J 
Since on the left-hand side of (1) we are summing binomial probabilities, and since 
the binomial distribution approaches the normal distribution as n increases, we have 
by the central-limit theorem that K = np + 8(n, d) (npq)'”, where for any given positive 
value A, a number 7, can be found such that if n > n, , then 8B, — A < B(n, A) < 
8, + A. We denote 8(n, d) (pq)'” by B. Then K = np + Bn™, where for n > ng, 
(8, — A) (pq)'” < B < (8B, + A) (pq)'”. 


Now 
. *) 


and using Stirling’s inequalities r’(2r)'”” exp (— r) < r! < r'(2 rr)'” exp (—r + 1/12r), 


1960 





we find 
/ (x) 
, | (2) 
K*(2rK)'” exp (—K+1/12K)(n—K)"*[2x(n—K)]'” exp {—n+K+1/12”—K)} ; 
n"(2rn)'’* exp (—n) 
Setting K = np + Bn”, simplifying the right-hand side of (2), and taking logs, we get 
—log (") 8 ————EEe 
\K i * " 12[npq + Bn'’"(q — p) — B*) (3) 


np + Bn’? + 3) log (p + Bn™'”) + (nq — Bn’? + 4) log (q — Bn™™”). 
Expanding log (p + Bn~’”*) and log (q — Bn~'””) around p and q respectively, and 


substituting into (3), we get 


—log (".) < n(p log p + q log gq) — n'”B log (q/p) + & logn, (4) 
K 


where © is a positive constant which does not depend on n. Converting the logarithms 


in (4) to the base 2, we find 


— log, ee) < n(p log, p + q log, q) — n'”[B log, (q/p) — Qn™” log, n]. 


K 
Then 
log, L < n(i + p log, p + q log, g) — n'”[B log, (q/p) — Qn7"” log, nj. 
If 
n >a, log, L < nC — n\[B\(pq)'” logs (q/p) — A(pq)'” logs (q/p) — Qn-"” log, n]. 


This proves the theorem, since A can be taken as close to zero as desired, and n™'”” 
log.n approaches zero as n increases. 

3. A theorem of Cramér. Before discussing the semi-continuous channel, we quote 
a theorem of Cramér, [1], which will be used later. 

Z, , Z, , +++ are independent, identically distributed chance variables, each with 











212 LIONEL WEISS (Vol. XVIII, No. 3 


probability density function v(z), expectation zero, and finite positive variance o”. 
There is a value A > 0 such that 


a co 


R(h) = | exp (hz)u(z) dz convergesfor |h| < A. 


m(h) denotes the first derivative of log R(h) with respect to h, and o{(h) denotes the 
second derivative of log R(h) with respect to h. We define A, as sup {h | R(h) converges}. 


From our assumption, A, is positive and may be infinite. M is defined as o~* lim,.4,-, m(h). 
M is positive and may be infinite. 
Theorem. For any value g in the open interval (0, M), the equation m(h) = og has a 


unique root h(g), which is positive. P(Z, + ++: +Z, > agn) = 
n~** {h(g)a,(h(g))(2m)'7 1? + Q(n, g)] exp {—n[h(g)m(h(g)) — log R(h(g))]} 


where for any values y; , y2 withO <7; < y2 < M, there ts a positive finite value B(y, , ¥2) 
with n | Q(n, g) | < Bly: , ye) for all g in the closed interval [y, , y2]. 


4. The semi-continuous symmetric channel. We assume that 


[ 8Goy) + Ay) Mog {2Kfo(y) + AIWF dy 
exists. Also, for any value y for which f,(y) = /,(y) = 0, we consider 2[f,(y) + f.(y)! 
0, 1, and 0 log‘0 is always to be considered as equal to zero, for any positive 


f(y) 


as zero for2 = 
k. The capacity C for the semi-continuous symmetric channel is 


a & 


| foly) logs {21fo(y) + A(T 'folw)} ay. 


For any given sequence (2% , , %,) of zeros and ones, W(x, , --- , 2) denotes a 


subset of (Y, , --- , Y,) space such that 

(a) P((Y,,-+: , Y,) isin W(a, , «++ , 2,)) > 1 — A, when the joint probability 
density function of Y, , «++ , Yn is fe,(ys) fes(Y2) ++ * fen(Ynd3 

(b) P((Y,,-::,Y,) isin W(z, , --- , x,)) is minimized when the joint probability 
density function of Y, ,--- , Y, is 


n 


[2(folys) + fi(ys))], 


t=] 
subject to (a). Then it is clear that the probability in (a) will be exactly 1 — ». We 
want to find the value of the probability in (b), which we denote by P*. 
By the familiar Neyman-Pearson lemma, the region W(z, , --- , 2,) is given by the 


set of points (Y, , --- , Y,) with 
I] {20fo(¥s) + ACY '¢(Y)} = k, 
i=1 


where k is a properly chosen constant. Alternatively, W(x, , --- , z,) is the set of points 


(Y,, 2 °* 5 Fg) With 


Y log {2ifo(¥) + f(YII%e(Y)} > log k. 
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By the definition of the symmetry of the channel, the distribution of 


p> log {2[f(Y.) + ACY #AY)} 
when the joint probability density function of Y, , --- , Y, is []?-, f.,(y:) does not depend 
on the sequence (x, , --- , 2,). Therefore the value of log k does not depend on the 
sequence (2, , -** , Xn). 
Also, the probability P* does not depend on the sequence (x, , --: , 2,). To show 
this, it suffices to show that the distribution of 2[f.(Y) + fi(Y)]"‘fo(Y) when Y has the 
probability density function 3(fo(y) + fi(y)) is the same as the distribution of 


2fo(Y) + f:(Y)])-'f:(Y) when Y has the probability density function 3(fo(y) + fi(y)). 
The ‘th moments of these expressions are respectively 


| 2" Toy) + ADIL) ay, [ 2) + AIO ay. 


But these integrals are respectively the (¢ — 1)st moment of the chance variable 
2f,(Y) + f:(Y)}* fo(Y) when Y has the probability density function f,(y), and the 
/ — 1)st moment of the chance variable 2[f,(Y) + f,(Y)]~’ f,(Y) when Y has the prob- 
ability density function f,(y), and therefore are equal. Since the moments determine 
the distribution of a bounded chance variable, the demonstration is completed. 
Therefore neither log k nor P* depends on the sequence (x, , --- , 2,), and it is no 
loss of generality to assume that x, = --- = x, = 0. We denote 


a 


| fo(y) log 2[foly) + fill 'foly)} dy by H, 


and 


[ fuyog {2Kfoly) + AI 'foy))? dy — H? by J?, 


J being taken as positive. Then, by the central-limit theorem, log k is equal to nH — 


n, d)n'’* J, where B(n, d) has the properties described in Sec. 2. 
We denote 


[ Mfoly) + fily)] log {2[fo(y) + hI 'foly)} dy by S, 
and 
[3 Ltf.(y) + f(y) og {2lfo(y) + fi(y)]'foly)})? dy — S? by o?,¢ 


being taken as agree To find P*, we use Cramér’s theorem of Sec. 3, with Z; = log 
Of (Ys) + f(YO] fo(VY)} — S, and ogn = nH — nS — B(n, \)n'” J, so that g = 
(H — S)/o — phe \) J/on'’’, First we must verify that g is in the open interval (0, 1/). 
It is easily shown that m(h) is equal to 


[_ BGol) + ADT FC) log (2) + HOT") ay 





i aca slain 


[| BG.@) + KTR ay 
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from which it follows that m(1) = H — S, and that m(h) is a strictly increasing function 
of h for positive h, except in the trivial case (which we exclude) where f,(y) = f(y) 
almost everywhere. It is easily verified that H — S is positive. This proves that g is 
in the open interval (0, 17), at least for large values of n. 

Our next task is to find h(g). We note that dm(h)/dh is continuous in a neighborhood 
of h = 1, and is equal to J° at h = 1. Therefore we have m(h) = m(1) + (h — 1)J? + 


e(h — 1), where (1/r) e(r) approaches zero as r approaches zero. The equation m(h(g)) = 
og becomes 


H — § + (h(g) — 1)J* + eh(g) —-1) =H -S— pen,r)In™, 


or 
h(g) = 1 — Bn, d)/(n'? J) + 6, 
where n’”” 6, approaches zero as 7 increases. 

Substituting this value of h(g) in Cramér’s theorem, we find that log P* = — n 
[H — B(n, ») Jn™’” + Q,], where n'”” Q, approaches zero as n increases. Denote by J’ 
the quantity that J becomes when logarithms to the base 2 are used in the definition 
instead of logarithms to the base e. Then log ,P* = — n[C — B(n, d) J’n''” + D,], 
where n'’’D,, approaches zero as 7 increases. 

If we have a code (u, , A;), «++ , (uz , Ax) of length Z and probability of error not 
greater than X, it follows from our discussion that P((Y, , --- , Y,) isin A;) > P*, when 
the joint probability density function of Y, , --- , Y, is [[?-; [4(fo(y:) + fily,))]. Since 
A,,-:+:,Az are disjoint, it follows that LP* < 1, or that log .L < n[C — B(n, d) J’n'” 
+ D,], where n'” D, approaches zero as n increases. Thus the stronger converse is 
proved for the semi-continuous symmetric channel. 

5. Acknowledgment. The author would like to thank Professor J. Wolfowitz for 
many helpful discussions of coding theory. 
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THE RAYLEIGH-RITZ METHOD FOR DISSIPATIVE 
OR GYROSCOPIC SYSTEMS* 


BY 
R. J. DUFFIN 
Carnegie Institute of Technology and Duke University 


1. Introduction. According to a well known principle of Raleigh, the natural 
frequencies of a conservative mechanical system can not be decreased when constraints 
are imposed on the system [1]. This principle is the basis for the Rayleigh-Ritz method 
for the approximate determination of the natural frequencies. 

Rayleigh treated vibration about static equilibrium, but he did not consider vibration 
about sleady motion. Thus, his analysis does not apply to vibration of machines with 
electric motors. The diffleulty is that gyroscopic forces have a different nature from 
forces arising from a potential. 

In a previous paper, the writer treated highly dissipative systems [2]. A relation was 
introduced termed ‘‘overdamping”’, which insures that the general solution of the equa- 
tions of motion is a sum of exponential decays. The decay values of overdamped systems 
are analogous to the natural frequencies of conservative systems. It was shown that 
Rayleigh’s principle applies to the decay values. 

In this paper, a simple transformation is introduced, which converts the equations 
of motion for a conservative system into the equations of motion for an overdamped 
system. Moreover, the addition of gyroscopic forces does not change this correspondence. 
By means of this device, it is seen that Rayleigh’s principle is valid for the natural 
frequencies of vibration of such conservative gyroscopic systems. 

2. Vibrations of conservative systems about static equilibrium. An example of 
the type of system considered by Rayleigh, is the triple pendulum, shown in Fig. 1. 





ai 


> 





Iq; 
Fic. 1. Triple pendulum. 


*Received June 1, 1959. Prepared under Contract DA-31-124-ORD-10, Office of Ordnance Research, 
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The potential energy is of the form c(q)/2, where 
c(q) = Crh — Cuts a C'ss% +r 2(C 1291 42 + (239293 + C3193.) - 


Of course, terms of higher order are neglected. The force components are obtained by 
differentiating the potential energy, so Newton’s equations for small free vibrations are 


Aig’ + >> Cig; = 0. 


Here, A; , A, , and A; , are the masses of the bobs, and a prime denotes differentiation 
with respect to the time, ¢. 
In matrix notation, the equations of motion may be written as 


Aq’ + Cq = 0. (1) 


In the general case of n degrees of freedom, g denotes a vector of nm components, and A 
and C are real symmetric matrices with positive definite quadratic forms a(q), and c(q). 
Otherwise, A and C are arbitrary. 

A normal mode of vibration is a solution of (1) of the form 


q= ee. (2) 


Here, u is a constant vector, and the constant w is termed a natural frequency. The 
equation which determines the normal modes of vibration is 


—w’ Au + Cu = 0. (3) 
A constraint is a linear homogeneous relation imposed on the coordinates such as 
q1 = 242 + q3 = 0. (4) 


For the triple pendulum, this constraint could be achieved physically by replacing the 
two strings connecting the bobs by a rigid rod. 

To study the effect of constraints on the spectrum of natural frequencies, Rayleigh [1], 
introduced the functional 

p (q) = = . (5) 

The stationary values of p are precisely the natural frequencies w; . The stationary values 
of p when q is restricted by (4) give the natural frequencies w/ of the constrained system. 
An additional constraint imposed leads to the frequencies w/’ of a doubly constrained 
system. The qualitative relationships are brought out in Fig. 2. Considering only the 
positive frequencies, we see that the constrained spectrum separates the original spectrum 
(This was proved independently by Routh, [3)]). 

3. Overdamped systems. The presence of dissipative forces can be accounted 
for by adding an appropriate term to (1). Then, the equation of motion becomes 


Aq” + Bq’ + Cq = 0. (6) 
Here, B is a symmetric matrix whose quadratic form b(q’) gives the rate of dissipation 
of energy [1]. 


Attention is now focused on systems where the dissipation is large. This would be 
the case, for example, if the triple pendulum were immersed in a very viscous fluid. 
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Fic. 2. Effect of constraints on the frequency spectrum. 


The precise mathematical condition to be assumed is that the quadratic forms satisfy 
the following additional restriction 
b?(v) — 4a(v)cv) > 0 (7) 
for an arbitrary real vector v ¥ 0. In a previous paper, this was termed the overdamping 
condition [2]. It was shown that such onerdamped systens have normal modes of decay. 
q = ue". (8) 
Here, u is a constant vector, and the real constant k may be termed a decay value. 
Substituting (8) in (6) gives 
k?Au + kBu + Cu = 0. (9) 
Forming the scalar product of (9) with u gives 
k’a + kb +c = 0, and 
k = [—b +(b® — 4ac)'”]/(2a) 
It was found that » of the normal modes require the positive choice of the radical in (10). 


These modes are termed primary, and their decay values are k, , k, , --- , k, . Likewise, 
n normal modes require the negative sign. These modes are termed secondary, and their 


decay values are h, , hz, +++ hn. 
The following variational principle may be seen to apply to overdamped systems. 


Let the primary functional be defined as 

pv) = [—bv) + {b°(e) — 4a(r)eo)}'*]/[2a@)]. (11) 
Then the primary decay values are the stationary values of p(v). Let the secondary 
functional be defined as 


sv) = [—bv”) — {b°@) — 4a)e(v)}'”"]/[2a(r)]. (12) 
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Likewise, the secondary decay values are the stationary values of s(v). 

A constraint on an overdamped system is seen to maintain condition (7). Thus the 
constrained system is also overdamped. The influence of constraint on the decay spectrum 
is illustrated in Fig. 3. A constraint is indicated by a prime. This diagram brings out 
the fact that the primary decay values of the constrained system separate the primary 
decay values of the original system. A similar separation relation holds for the secondary 
spectra. 

It is seen from Fig. 3 that the following inequalities hold 


h, < sv) <h; < ki < pve) < ks < 0. (13) 


primary 


secondary 





Fia, 3. Effect of constraints on the decay spectrum. 


Here, v is an arbitrary vector; it results from a double constraint. Clearly, (13) would be 
useful for estimating the decay values of the overdamped triple pendulum. 

These considerations show that the analogue of Rayleigh’s principle holds for an 
overdamped system, provided the primary spectrum and the secondary spectrum are 
considered separately. Of course, even in the classical case, shown in Fig. 2, it is necessary 
to distinguish between the positive frequency spectrum and the negative frequency 
spectrum. This analogy is not accidental. Thus, by a trivial relaxation in the definitions, 
it is now to be shown that a conservative system is a special case of an overdamped 
system. 

A relative overdamped system is defined by matrices A, B, and C, such that: 


(i) They are symmetric. 

(ii) The overdamping condition (7) holds. 

(iii) A is positive definite. 
With A, B, and C having these properties, let k = ky + m, in Eq. (9). Then (9) becomes 
(27Au + k,Bou + Cou = 0. Here By = 2mA + B, and C, = m’A + mB + C. It may be 
verified that 


be — 4ac, = b° — 4ac > 0. 
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Moreover, if m is a sufficiently large positive constant, both B, and C, are positive 
definite. It follows that A, Bo , and Cy define an ‘‘absolute” overdamped system. Thus, 
the spectrum of a relative overdamped system is simply the spectrum of an absolute 
overdamped system shifted upward. 

The normal modes of vibration of a conservative system are defined by Eq. (3). 
But, clearly, the matrices A, 0, —C, correspond to a relative overdamped system. This 
proves that the conservative system is a special case. 

4. Vibrations of conservative gyroscopic systems. Now of concern are mechanical 
systems whose small motions are governed by an equation of the form 


Aq’ + Bd + Cq = 0. 

The matrices are to satisfy the conditions: 

(i) A and C are symmetric positive definite. 

(ii) B is antisymmetric. 
We term such systems regular. It is seen that regular systems are conservative because 
the quadratic form b(q’) vanishes identically. 

\ simple example of a regular system is the gyroscopic pendulum shown in Fig. 4. 
A bob of mass m is spinning so as to have angular momentum J, directed along the rigid 
supporting rod of length h. The position of equilibrium is taken as origin of x, y coordinates 


in a horizontal plane. 


Fic. 4. Gyroscopic pendulum. 


Taking moments about the point of support, leads to the following equations of 


motion. 


mhx’’ + Lh y’ + mgz 


mhy’’ — Lh7'x’ + mgy 
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It is easy to vertify that the normal modes of vibration are rotations in a circle with 


— * sie 2 ee 
o = mht = | (=Z43) + 2| ‘ 


The positive sign corresponds to counterclockwise rotation. A complete analysis of such 
pendulums is given by A. G. Webster [4]. General gyroscopic systems are treated by 
E. T. Whittaker [5]. 

teturning to the general case, it is seen that the equation for the normal modes of 


frequencies 


vibration of a regular system is 
—w’ Az + iwBz + Cz = 0. (14) 


Here, the vector z may be complex, but it can be written as z = u + w, where u and v 
are real vectors. This expression for z is substituted in (14), and the real and imaginary 
terms are separated. The two equations resulting, may be written as the following single 
equation, using compound matrices. 


2 AO OB CO U _ , 
[=(04) + (0) ~ (6c) (2) = @ as 


The matrices in (15) are symmetric. The overdamping condition is satisfied because 
4(u-Bv)? + 4[a(u) + a(v)][c(u) + cv)] > 0. 


Moreover, the first matrix in (15) is positive definite. Thus, (15) may be regarded as 
the equation for the normal modes of an overdamped system with 2n degrees of freedom. 

If w and (u, v) satisfy (15), then w and (v, —u) also satisfy (15). Thus, each decay 
value of (15) is “double”. Otherwise, the decay spectrum of (15) is identical with the 
frequency spectrum of (14). 

A constraint on (14), corresponds to two constraints on (15), one on u, and one on v. 
Or course, the constrained spectrum of (15) is double valued. Then it is seen that in- 
equalities of the type (13) lead to the following separation property. The double values 
of the doubly constrained spectrum separate the double values of the original spectrum. 

By the correspondence between the spectra of (14), and (15), it follows that the 
separation property holds for the natural frequencies of (14). Moreover, if w is a natural 
frequency, it is seen from (14), that —w is also. Thus, the frequency spectrum is exactly 
as shown in Fig. 2. This proves Rayleigh’s principle for regular systems. 
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Fig. 5. Network analogue of the gyroscopic pendulum. 


5. Non-reciprocal electrical networks. The classical analogy between mechanical 
systems and electrical networks makes force correspond to voltage, and velocity cor- 
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respond to current. The matrices A, B, and C, correspond to matrices of inductance, 
resistance, and (capacitance)~*. If A, B, and C, are symmetric, the network obeys the 
reciprocal law of Rayleigh [1]. 

A recent development in the network art has been networks which are the analogue 
of conservative gyroscopic systems. To aid in the analysis and synthesis of such non- 
reciprocal networks, Tellegen [6], introduced an ideal network element termed a gyrator. 
This concept is brought out in Fig. 5, which gives a network analogue of the gyroscopic 
pendulum. 

It is seen that the gyrator serves as a representation of the antisymmetric B matrix. 


REFERENCES 


1. Rayleigh, Some general theorems relating to vibrations: Sec. I. The stationary condition, Sec. II. The 
dissipation function, Sec. III. A law of reciprocal character, Proc. London Math. Soc. 4, 357-368 (1873); 
Scientific papers, 1, pp. 170-184 

. R. J. Duffin, A minimaz theory for overdamped networks, J. Ratl. Mechanics and Anal. 4, 221-233 
(1955 

3. E. J. Routh, Dynamics of rigid bodies, vol. 2, 6th ed., Macmillan, London, 1905, p. 57 
A. G. Webster, Dynamics, 2nd ed., Stechert, New York, 1922, pp. 288-296 

. E. T. Whittaker, Analytical dynamics, 4th ed., Macmillan, 1937, pp. 191-195 

. B. D. H. Tellegen, The gyrator—a new electric network element, Phillips Research Reports 3, 81-101 

1948 


or 


a 











BOOK REVIEWS [Vol. XVIII, No. 3 


222 


BOOK REVIEWS 


Mathematical programming and electrical networks. By Jack B. Dennis. Technology 
Press of The Mass. Institute of Technology, John Wiley & Sons, Inc., New York, 
Chapman & Hall, Ltd., London, 1959. vi + 186 pp. $4.50. 


This well-written and well-organized little volume, the third in the MIT Research Monograph 
Series, presents the results of the novel and recent researches of the author in establishing an equivalence 
between mathematical programming, particular linear and quadratic, and electrical networks composed 
of voltage and current sources, ideal diodes, resistors, and transformers. 

The central idea is based on the extremum principle that in a circuit of current and voltage sources, 
resistors, and transformers the currents and voltages in the various branches of the network are such 
as to minimize the heating power losses. This is essentially a minimization of a quadratic form in the 
currents when resistors are present and a linear form when no resistors are present. Additional linear 
constraints are provided by the Kirchoff node laws. The introduction of diodes into the circuit imposes 
a condition of non-negativity on the currents. Thus, with resistors one has the electrical equivalent 
of a quadratic programming problem and without them, a linear programming problem. 

The first two chapters and the Appendix, which occupies the last third of the book, provide some 
general mathematical notions for mathematical programming, local and global minima, concavity 
and convexity, Lagrange multipliers, the Legendre transformation, duality, and a number of pertinent 
theorems. In view of the greater advantages, in speed, economy, accuracy, and facility, of the use of 
large scale digital computers in solving these problems, as the author has pointed out, it is unlikely 
that anyone would use these techniques to solve a practical problem. Nevertheless, this research provides 
essentially a new outlook and insight into mathematical programming which has heretofore relied 
almost completely on economic, scheduling, and allocation examples as the vehicles for the introduction 
to mathematical programming. The minimum employment of the economic language and examples 
so common to other introductions to linear programming should be to most eng’neers’ liking. 

Chapter Three describes the properties of the various electrical devices of the net works, the network 
laws, and the equivalences of the various network problems. This possibly should appear a little earlier 
so that the non-engineer can see what the author is driving at. Chapter Four presents a diode-source al- 
gorithm for solving electrical models of certain network flow problems. For the maximum flow problem 
this is similar to the Ford-Fulkerson algorithm for transportation problems. The central notion here, 
useful throughout the remainder of the book, is that of a breakpoint curve, a curve of possible voltage- 
current combinations for a diode. The algorithm, elaborated further in Chapters Five and Six, consists 
in following this curve through in a manner analogous to the consecutive steps of the simplex method 
of Dantzig. Certain modifications to the algorithm and several other electrical algorithms called ‘‘valve’’ 
and “by-pass” algorithms are also presented and applied to different problems. 

The final chapter proposes an iterative method for solving general equality- and inequality-con- 
strained minimization problems based on a quadratic program. This method does not appear to differ 
widely from some methods in current practice. Errors are very few and relatively minor in consequences. 


The bibliography is adequate. 
M. L. Juncosa 


An introduction to statistical mechanics. By J. S. R. Chisholm and A. H. de Borde. 
Pergamon Press, Inc., New York, London, Paris, Los Angeles, 1958. ix + 160 pp. 
$6.00. 

This book is designed as a text suitable for an Honours course in British Universities. The subject 
matter is limited to equilibrium phenomena and the applications are those which can be found in most 
books on statistical mechanics written during the last twenty-five years or so. The main virtues of the 
book are that it is short and clearly written. The scope, however, is very limited and there is no indi- 
cation of the modern trends of the subject. The author claims to have a simplified version of the Darwin- 


(Continued on p. 234) 
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FURTHER PROPERTIES OF CERTAIN CLASSES 
OF TRANSFER FUNCTIONS* 


BY 
ARMEN H. ZEMANIAN 


College of Engineering, New York University 


Abstract. Some previously defined classes of rational transfer functions are extended 
to transfer functions that need not be rational. Several additional properties of these 
types of functions are then developed. Finally, bounds on certain derivatives of the unit 
impulse responses corresponding to such functions are shown to exist. 

1. Let Z(s) denote the transfer function of a system and let W(t) be the response 
of the system to a unit impulse applied at time ¢ = 0. Then Z(s) and W(?) are related 


by the Laplace transform. 


Zs) = [ ” W(t) exp (—st) dt. (1) 


0 


It will be assumed henceforth that W(¢) and its derivatives are integrable in any finite 
interval. The transfer function Z(s) is a function of the complex variable, s = ¢ + jw. 
Since W(¢) is a real function for physical systems, the complex singularities of Z(s) 
occur in complex conjugate pairs. 

In a recent series of papers [1, 2, 3] a number of properties were established for 
certain classes of rational transfer functions. One of the objects of this paper is to remove 
the restrictions that Z(s) be rational. This will require a more general form for the defi- 
nitions of the classes of transfer functions. Let Z(s) be analytic for ¢ > 0 and let it be 
asymptotic to K/s* as s approaches infinity where K is a positive constant and k is 
a positive integer. In the special case where Z(s) is rational, we have 








— 8” + 4,8" +++ Ms 
Zs) = K s" + bs" + +++ + dy ” 
k=m-—n>O0O. 8) 


As noted above, W(¢) assumes only real values. Consequently, all the coefficients in (2) 


must be real. 
Returning to the general case where Z(s) need not be rational, let Z,(s) denote the 


following successive integrations of Z(s). 
Zo = / ds... / di ses / Z(&) di « (4) 


Only the principal branch of this multivalued function will be needed. In fact, it will be 
sufficient to assume that the arbitrary paths of integration in (4) never enter into the 
region defined by « < o, < 0 where og, is the largest real part among the singularities 
of Z(s). Under this restriction, (4) can be considered a single-valued function. 





*Received June 19, 1959. 
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Also, let the real and imaginary parts of some subsequently needed quantities be 
denoted by 


s=ao-+ p, (5) 
8, = og + hy ; (6) 
Z(jo) = R@) + jl), (7) 
Z (jo) = Rw) + j1,). (8) 

The following relations may then be obtained from (4). For q odd, 
R,@) = (-1)*?” [ dug i diy-2 *** [ I (wo) duro (9) 

and, for g even, - aa . 

Ril) = (—1) [ d..1 [ Py ee [ Rw.) hie (10) 


R,(w) and J,(w) are even and odd functions of w, respectively [2; lemma 2]. 
We may now state the following definition. 
Definition 1. A function Z(s) will be called a class k function if the following conditions 
hold. 
(a) Z(s) zs analytic for o > 0. 
(b) Z(s) ~ K/s° as s — © where K is a positive constant and k is a positive integer. 
(c) Z,-1(s) ts a positive real function in the half plane, o > 0. 


A function is said to be positive real if its real part is nonnegative for o > 0 and 
if it assumes only real values along the real axis [4]. 
A certain subclass of each class k will be required in the following discussion. The 


definition of this subclass is again a generalization of the previously given one [2] in 
that the condition that Z(s) be rational is dropped. Furthermore, R(w) and d//dw are 
now allowed to assume the value of zero at w = 0. 


Definition 2. A function Z(s) will be called a subclass k function if conditions (a) and 


- 


(b) in Definition 1 hold and if the following condition holds. 
(c’) For k odd, R(w) has k — 1 changes of sign for — ~- <w< @ and R(w) ts positive 
in the neighborhood of w = 0; for k even, I(w) has k — 1 changes of 
sign for — © <w< o and dI/dw is negative in the neighborhood of w = 0. 


It has been shown formerly [2; theorem 2] that all subclass k functions are class k 
functions. Even with these more general definitions, the previously given proof applies 
in precisely the same way. Furthermore, all the theorems in Part II of Ref. [2] and in 
Ref. [3] continue to hold since their proofs did not make use of the rationality of Z(s). 
(Theorem 6 of Ref. [3] must be modified slightly if it is to apply in this case.) 

One particular result [2; theorem 9] that the unit impulse response corresponding 
to any class k function must satisfy is 

Ki 


I\wols &—D! 


DI (11) 


The principal objective of this paper is to show that the following bounds exist on a 
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number of the derivatives of those W(¢) corresponding to subclass k functions. In par- 
ticular, fork > 2 


rw) a 
[werols i. a (12) 


where the integer u is restricted by 0 < u < (k — 1)/2 if k is odd and by 0 < » < k/2 
if k is even. In developing this result, some new properties of the subclass k functions 
will be obtained. 

2. The following lemmas will be needed. 

Lemma 1. Let F(s) satisfy conditions (a) and (b) of Definition 1. Then, for k odd, 
R(w) must have at least k — 1 changes of signin — ~ < w < @ and, for k even, I(w) 
must have at least k — 1 changes of signin — © Cw < @, 

The proof of this lemma has been given previously [2; lemma 3]. It applies in the 
case where F(s) is not rational since the rational property was not invoked anywhere 
in its proof. 

Lemma 2. Let k be even and let Z(s) be a subclass k function. Then sZ(s) is a subclass 
(k — 1) function. 

Proof. By condition (c’), I(w) has k — 1 changes of sign in — © < w < o. Since 
I(w) is an odd function [2; lemma 2], one of these sign changes occurs at w = 0. Now 
the real part of sZ(s) equals — w/(w) and, consequently, it has (k — 2) sign changes in 
— «© <w< o., Since dJ/dw is negative in the neighborhood of w = 0, — wl(w) is 
positive in the same neighborhood. Thus, — w/J(w) fulfills condition (c’). Conditions 
(a) and (b) are also fulfilled by sZ(s) so that the lemma is established. 


Theorem 1. Let k > 3 and let Z(s) be a subclass k function. Then, 


© 


G(s) = [ 8Z (So) dso (13) 


is a subclass (k — 2) function. 


Proof. The theorem will first be established in the case where k is odd. Denoting 
the real part of G(jw) by Rg(w), we have 


a i] " ile (14) 


Since sZ(s) is analytic for ¢ > 0 and since Z(s) = 0 (1/s") as s > @, (14) assumes the 
value of zero when the lower limit is set equal to — o. Hence, 


R,(w) = [ wWR(wo) dw . (15) 
Integrating (15) by parts, the value of Rg(w) at w = 0 may be expressed as 
R,(0) =_-_ [ dw, il R(wo) dw ° 


Because Z(s) is a subclass k function, f°. R(wo) dw will have k — 2 changes of sign in 
— © <w< © one of which is-at the origin. Hence, 
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/ hi [ Rw) devo (16) 


will have k — 3 sign changes in — © < w < o., It follows that, since R(w) is positive 
in the neighborhood of w = 0, (16) is negative at w = 0 so that Rg(0) is positive. 

Now consider the case where k = 4v — 1 (v = 1, 2, 3, ---). The quantity wR(w) and, 
consequently, Rg(w) become ultimately positive as w decreases indefinitely toward 
— o, Hence, for Rg(0) to be positive, Rg(w) can have only an even number of sign 
changes in — © < w < 0. However, wR(w) has (k — 1)/2 sign changes in — ~ <w < 0 
and this is an odd number. Since integrating according to (15) can never increase the 


number of sign changes in — © < w < 0, Re(w) must have less than (k — 1)/2 sign 
changes in — ~ < w < 0. Moreover, by lemma 1, Rg(w) has at least (k — 3)/2 sign 
changes in — ~ < w < Q. Thus, Rg(w) has exactly (k — 3)/2 sign changes in 
—-2x7<w< 0. 


This result coupled with the facts that Rg(w) is even and R¢(0) is positive shows 
that G(s) satisfies the requirements of condition (c’). Conditions (a) and (b) are also 
fulfilled so that G(s) is a subclass (kK — 2) function. 


A similar argument may be applied in the case where k = 4v + 1 (v = 1, 2, 3, ---). 
Now let k be even. As before, it can be shown that J¢(w) is given by 
Toles) = [col (en) des . (17) 


This quantity is continuous and an odd function of w. Hence, Jg(0) = 0. Also, its de- 
rivative is negative in the vicinity of w = 0 since J(w) has this property and J(0) = 0. 
Since Z(s) is a subclass k function, w/(w) has (k — 2)/2 changes of sign within the 
interval — ~ < w < 0. Let w; be a point where Jg(w) changes sign. Since J¢(w) is related 
to wI(w) according to (17), Z¢(w) must have a smaller number of changes of sign in 
— © <w < w; than w/(w) does. Furthermore, the point w = 0 is a point where J,(w) 
changes sign so that J,(w) has no more than (k — 4)/2 changes of signin — ~ <w < 0. 
Hence, invoking Lemma 1, Jg(w) has exactly this number of sign changes in — ~ < w < 0. 
It follows that G(s) is a subclass (k — 2) function. This completes the proof. 


i 


y, 


Theorem 2. Let Z(s) be a subclass k function where'k > 2. Also, let u be any integer in 
the range 0 < uw < (k — 1)/2 tf k ts odd and let u be any integer in the range 0 < uw < k/2 
if k is even. Then s*Z(s) is a class (k — pu) function. 


Proof. If a function F(s) is a class k function, then the quantity (— 1)*d'F/ds" is a 
class (k + h) function. This result is an immediate consequence of the definition of a 


class k function. 
Let k be odd (k = 2 + 1; v» = 1, 2, 3, ---). To establish the conclusion, it will be 


shown that the function, 


[ as,-s | ee / SZ (8) d8o , (18) 


a 


is a class (k — 2u) function where 0 < »w < vy = (k — 1)/2. 
Consider the innermost integral of (18). Integrating by parts repeatedly we may 


write 
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[ s52(60) ds, = s** f 52(6s) ds + (u— 1) [st Pads, | 50%(6x) do 


s 1 


[ &Z(8>) dS» 


vs 


+ (u — Ist? | ®, de, i] $Z(8) a8 


1 


Te / si ds, / o, de, / &Z(8>) a8» « 


Continuing this process of integrating by parts the last term in this sum, the following 
result may be obtained. For yu odd, 


34Z(s.\ ds, = apt [ SoZ (So) ds» - cee 


(19) 
+ (u — 1)(u — 3) «++ 4-2 [ Sona @o-un.*** [ 8 Z(8) d8o 
and, for » even, ° 
[ s:Z%s,) ds =" [ 850(s,) diy + ++ 
= sid (20) 
+ (u — 1)(u — 3)-:-- aa-f ds,/2 [ S,o-1 A8yye-1 ° °° [ *sotis. dS» . 
‘ outs - 


Repeatedly integrating in the above manner, (18) may be rearranged_into the following 


finite sum. 


| Us,—3 | ds,-2 ote [ 8}Z (So) dso 
| Op-l ds,,. 1 [ 8,-2 ds,-» sia [ 89Z (Sp) ds» 
> 7 8yu-1 J 8, (21) 
+ A; ds, [ ds,_ [ 8,-2 8,2 **° [ 8.Z (8) d8o 


as A» “dees [ as, [ ds, [ ds,-2 | §,-3 ds,-3 Pee / &Z (80) dS» 


@ But. J 8y—1 


In this expression, all the A; are positive integers. 

Now Theorem 1 may be applied repeatedly to each term on the right hand side of 
(21) to show that each term is a subclass (k — 2u) function. That is, working from the 
innermost integral outward, Theorem 1 indicates that, when the integrand possesses 
the factor s, each integration yields a function of some subclass and the order of the 
resulting subclass is reduced by two from the order of the subclass of the initial function. 
Moreover, when the integrand does not have the factor s, each integration produces a 
subclass function whose order is reduced by one. Consequently, each term on the right 
hand side of (21) is a subclass (k — 2u) function and, hence, a class (k — 2) function. 

Furthermore, if F(s) and G(s) are class k functions and if A and B are positive 
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numbers, then it follows from the definition of the class k functions that AF(s) + BG(s) 
is also a class k function. Therefore, (21) is a class (k — 2y) function and, by the remarks 
in the first paragraph of this proof, s"Z(s) is a class (k — yu) function. 

The same result may be established in the case when k is even. (kK = 2p; 
vy = 2,3, 4, --- . Lemma 2 states the theorem in the case when k = 2.) The same pro- 
cedure as used before may be applied to the expression, 


ao a2 


dsj. | dss +* | StZ(s>) do , (22) 
J 8yu-—32 81 


Js 


to show that it is a class (k — 2u + 1) function where 0 < pw < v = K/2. In this case 
(22) is rearranged into the following finite sum. 


- rp a 


s | 8-9 ds,,-» | a8 ds, Te | 8 Z (So) ds, 


nx x + 


+ B, | ds,-y [ ae | SoZ (So) dS 


In this expression, the B; are all positive integers. Moreover, each term is a subclass 
I g , 
(k — 2u + 1) function by the same argument as before. Hence, (22) is in class 
(k — 2u + 1). Differentiating it 1 — 1 times and multiplying by (— 1)*~* will produce 
, 5S J 
a class (k — yu) function. This completes the proof. 
3. The principal conclusion of this paper may now be stated. 


9 


Theorem 3. Under the hypothesis of Theorem 2, 


"ii 
| w(t) | < K —————- 
, OIS (k — p— 1)! 

Proof. The transfer function Z(s) is asymptotic to K/s* (k > 2) as s —> ~. Hence, 
by the initial value theorem [5; theorem 15, p. 267], the corresponding unit impulse 
response W(t) and its first k — 2 derivatives all have an initial value of zero. Thus, for 
the stated ranges of u, s*Z(s) is the Laplace transform of W“ (¢) [5; p. 129]. 

Finally, it has been shown [2; theorem 9], that the unit impulse response correspond- 
ing to a class k function is bounded according to (11). Consequently, the conclusion 
follows immediately from Theorem 2. 
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ON THE NON-EXISTENCE OF TRANSONIC PERTURBATIONS* 


BY 
J. D. NEETHLING 
National Physical Research Laboratory, Pretoria, South Africa 


M. Schafer [2] has developed a practical method to construct plane transonic potential 
flows directly by introducing a mesh of derived characteristics so defined that the mesh 
is real in both the supersonic and subsonic regions and joins continuously on the sonic 
line. Since the general solution of the equations for derived characteristics leads to great 
mathematical complications, Schafer considers a simple solution which leads to the 
well-known Chaplygin solutions of the gasdynamic equations. These are a one parameter 
family of solutions which by suitable choice of the parameter, can be made to describe 
flows past profiles of practical interest. 

Since the equations are linear in the derived characteristic variables which, in the 
case of the simple solution mentioned above, coincide with the polar hodograph variables, 
it is possible to find new solutions by superimposing solutions with different values of 
the parameter. In a second paper Schifer [3] discusses the perturbation of profiles 
corresponding to special values of the parameter by superimposing solutions with 
positive integral values of the parameter on these specially chosen solutions which 
may be regarded as basic flows. The superimposed solutions are multiplied by coefficients 
which can be arbitrarily chosen to satisfy perturbed boundary conditions. Schafer 
studies perturbations of the boundary conditions which are given by analytic (ex- 
ponential) functions but which are practically localized in the supersonic region. 

C. S. Morawetz [1] on the other hand, shewed that no perturbations involving only 
a finite segment of the supersonic boundary of a plane potential flow exist. If the per- 
turbations studied by Schiifer could be regarded as numerical approximations of strictly 
localized perturbations, we would indeed have a family of flows depending—not 
necessarily smoothly'—on a parameter and with the same velocity at infinity namely, 
zero, which differ only along a finite segment of the supersonic boundary. Moreover, the 
theorems of Morawetz cannot be applied directly to these flows since in general a direc- 
tion singularity occurs at infinity and hence infinity cannot be represented by a point 
in the hodograph plane considered by Morawetz. 

irstly, we remark that if the superposition is truncated, the perturbations must have 
an analytic character since each term in the summation is an analytic function of its 
variables. This means that singularities of curvature are impossible and the perturbations 
cannot be restricted to a finite segment of the supersonic boundary. We shew that the 
coefficients in the infinite series cannot be so chosen that a singularity is introduced 
into the boundary. 

Using the simple solution of the equations for derived characteristics mentioned 
above, Schafer [2] obtained the one parameter family of solutions of the gasdynamic 
equations in terms of certain metric functions as dependent variables and the polar 


*Received October 16, 1958. Revised July 7, 1959. 
1A set of flows given by infinite series of the solutions with arbitrary coefficients may be constructed 
by summing over the parameter. 
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hodograph variables v and # as independent variables, known as the Chaplygin solu- 


tions. In Schifer’s notation the solutions are [2]: 





is K 77?) ceed or (1a) 
p Leos kd 
St ie 1 [ko 'H,(v") + Qv** Hi (y")] Jos kd sl (1b) 
P Lsin ko 
k+l kal 
T— 2% = E 1 Hi sin (k + 1)8 + a (kH, + v°H}) sin (k — vo |, (2a) 
pLk + iia 
l “ 
oo 1] =v Hj, cos (k + 1) + oe 7 H, + v°H}) cos (k — V9 |, (2b) 
pLk 
, , » |coskd or 
y =vH,v) 4 : (3) 
sin kd 
where 
H, = F(a, — 8; y; v°), the hypergeometric function, 
and 
1H ,(w) 
Hap) = Cale), 
Hi(u) dw 


= and H are the metric functions giving the relation between the cartesian coordinates 
x and y in the field of flow, and the derived characteristics which in the case considered 
coincide with the polar hodograph variables v and # (v is the absolute value of the velocity 
vector made dimensionless with the velocity at zero pressure). y is the stream function, 
p the density made dimensionless with the stagnation density and «x the adiabatic 
gas constant. 

We shall require an asymptotic expression for the hypergeometric function, H,(v’), 
for large positive values of k. We write 

(ty) Taet+aTi~ sv" 


= - - v ( 
H, p> Iiy+2 T@ I-8) 2!’ (4) 


1 ra i ia, 
i * oe * , 
I J ie, 
1°¢=)’ 


andy = k+1.a> 0and8 > 0 for positive k. We have for large k 


My) Ta + 2) 
4 ; a oth 


Try + ae I'(a) 





where 


mn 
ae 
ero 
% 
Ht 


: f, 
\ 


] 
2 


| 
| 
a 
a 
| |+ 


and 


HW) ~ EB 5 salt (+ (EY) = [1-3 +E) YI: 
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By suitable choices of k Schifer makes the solutions (1) represent basic flows of 
practical interest. Perturbations are then applied to these basic flows by putting 


=, 8) = E*@, 8) + : > o'H.(’)(a, cos kd — b, sin kd), (6a) 
kok, 

Hiv, 8) = H*(v, 8) + : > (ko "Hi, + 2v**'HD(a, sin kd + b, cos kd). (6b) 
k=k, 


The symbols with stars denote quantities in the undisturbed flow and plain symbols 
quantities in the disturbed flow. Since the equations for Z and H are linear, the perturbed 
metric functions as defined by (6) are again solutions of these equations. It is clearly 
more convenient to satisfy the perturbed boundary conditions in the physical plane and 
the expressions for x and y can be supplemented accordingly since they are linear combi- 
nations of = and H. We shall now shew that the perturbation terms in (6) must be 
analytic functions of v and 3. We state the proposition in the following form: suppose 
the perturbation terms in (6) converge everywhere in the field of flow except possibly 
in a finite number of isolated points in the supersonic region. Then these terms converge 
absolutely and uniformly in all domains and they have term by term derivatives of all 
orders with respect to both variables. 
Consider the function 


vo *H,(v’). (7) 
This function tends to zero with increasing k, for in (5) v < 1 and 8 is asymptotically 
of the form 


bk + Ok), b>. 


The asymptotic expression for (7) is 
with 


V increases with v in the interval 


0<Yr< o{(‘ rs 1) iF + (s+) } 
kK — 1 kK — 1 


In the latter point V has a maximum whereafter it decreases to the value zero for 


1/2) -1 
vais (tt) =v. 
k—1 


We can therefore say that there is a limiting point of maxima of (7) at 


— g+1 1/2 < + ny _ 3 
v= 4(+1 1) E + (s =a vj; 


with respect to positive k. Since the treatment of the perturbation term in (6) is the 
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same whether we regard it as a sine series, a cosine series or a mixed series, we shall 
simplify by putting b, = 0 for all k. We now put 


a, = f(k)A* 


and derive restrictions on A and f such that the conditions of our proposition are satisfied. 
A is simply a normalizing factor for the expression (8). Let A be chosen such that 


1 p : 1 1/2 6 
Av,| 1 — 1 al at (<+2) )| = Es V4 cad Vv; . 


With v, in the range 0 < v4 < v, , the perturbation term in (6) would diverge as cosine 
series for all values of v in the interval v4 < v < v, orv, < v < vy according as v4 < v, 
or v; < v4". If the conditions of our proposition are still to be satisfied we must therefore 
have v4 = v; . On f(k) we must then obviously impose the restriction 


f(k) = ofA‘). 


The perturbation term can have at most one singular point with respect to v viz., the 
point » = v, . The order of this singularity is determined by the form of f(k). It follows 
that the perturbation term is a uniformly and absolutely convergent cosine series with 
v in any subinterval of the interval 


. 1/2) -1 
0<v’< of ~ (c+) = Us 
— 


which excludes v = v, . We exclude the range 
i 1/2} -1 
of 1 + (+4) } <v’ <1 
\ > — 1 


since in this range (8) becomes complex. We shall see that the mapping of the hodograph 
plane on the physical plane ceases to be unique at the point v = v,. 

We moreover know that each term of the cosine series is a uniformly and absolutely 
convergent power series of v in the range considered. Hence the perturbation term can 
be rearranged as a power series in v which has term by term derivatives of al! orders with 
respect to v in all subintervals which do not include the point v = 2, . 

Differentiating the perturbation term formally n times with respect to 3, we obtain 


cos kd 


>(—)'k"f(k) A'v' "A, or, with ¢ = [3n] + 3. 


(sin kd 
But in this expression we still have 
(—)'k"f(k) = 0(A*) since f(k) = o(A’). 


Hence the perturbation term also has term by term derivatives of all orders with respect 
to 3 unless v = », . Note that if the perturbation term is singular at this value, then 
there will in general be a singular line in the supersonic region viz., the isotach v = v, . 
We now shew that with the above restrictions on a, our proposition also holds for 
(6b). The asymptotic expression for the coefficients of the perturbation term in (6b) is 





1For we can find a kp for every such v such that A*v*H;, ~ 7 for all k > ko with n > 1. 
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+1 ‘ 1/m + 1\'7| (ek +1\% : 

et —pe(1+ GE) pix fie iF} ” 

The first factor is the same as for (6b) while the second factor has a zero in the point 

» = v, where V(v) has a maximum. (9) therefore increases steadily until it approaches 

v = v, where it drops off sharply to negative values with the same amplitude in the 

range v; < v < v,. The coefficients of (6b) are thus seen to have a limiting point of 

maxima in v = v, with respect to k. This limiting point is clearly never reached and 
away from this point the coefficients of (6b) become asymptotically 


{aofs - bof +H)" TF 


which is the same as the equivalent expression for (6a). Hence there can be no singularity 
in the perturbation term in (6b). 

If we require furthermore that singularities only be propagated along characteristics, 
we must also exclude the singularity in (6a) since an isotach cannot coincide with a 
characteristic. 

In order to exclude the range v. < v < 1 we consider the mapping of the hodograph 
plane on the physical plane. The Jacobian of this mapping is (see [3]) 


(x,y) _ 1 (* ny? ) 
av, 9) ov ra +e} (10) 


c 
By comparing (10) with (8) and (9) we see that the expression 


[500 + GY") 


is asymptotically a factor of the Jacobian so that the mapping will cease to be unique 
when v = v, . Moreover a limiting line will generally have appeared in the basic flow 
before this velocity is reached—this is e.g., the case in the basic flow with k = — $ 





considered by Schiéifer. 
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BOOK REVIEWS 


(Continued from p. 222 


Fowler method for obtaining the canonical distribution. The details of this, which are given in an ap- 
pendix, are sound but the outline of it given in the text is deceptive. In essence the book is a condensed 


version of the books by Fowler. 
G. H. NEWELL 


Theory of value, an axiomatic analysis of economic equilibrium. By Gerard Debreu. 

John Wiley & Sons, Inc., New York, and Chapman & Hall, Ltd., London, 1959. 

xli + 114 pp. $4.75. 

This book is a landmark in the development of mathematical economics. Coordinating the work 
that has gone into the central body of economic theory since Cournot and recasting it in the mold of 
contemporary mathematics, this slender volume gives the definitive treatment of the Theory of Value 
constructed rigorously from first axioms. 

Chapter 1 assembles all the mathematical concepts and propositions to be used in the economic 
chapters. Evidently a great deal of thought has gone into such organizational matters as notation and 
layout, and this pays off handsomely in clarity and readability throughout the book. The departure 
from traditional calculus-bound methods is apparent in this box of tools: set and point set theory, 
ordering, metric spaces, and topology including fixed point theory. Chapter 2 clarifies all economic 
concepts employed and will be particularly useful to the mathematical reader as a grammer of the 
economist’s language. Economists on their part will be delighted at the crystal clarity with which the 
basic concepts and underlying assumptions are here exposed. The analysis gets under way in chapter 
3 with a discussion of the theory of production. The axioms on technology are presented and discussed 
in order of decreasing plausibility. These include the technology assumed in Linear Programming 
models as a special case. The implications of profit maximization are studied and ‘‘comparative statics’’ 
is analyzed. The place of derivatives and of Lagrangian Multipliers is taken here by supporting hyper- 
planes to convex sets. The 11 pages of this analysis are the most compact complete treatment of the 
general theory of production in existence. 

In Chapter 3 the theory of consumption is constructed analogously on the basis of complete, transi- 
tive, continuous preference orderings on closed, connected consumption sets. The existence of (ordinal) 
utility functions is proved and the theory of consumers’ choice subject to wealth constraints is de- 
veloped. The high point of the book is reached in Chapter 5 where the existence of an equilibrium for 
a “private ownership economy”’ characterized essentially by non-increasing returns in production and 
insatiable consumption subject to non-increasing marginal utility (these are old-fashioned economic 
terms not used here) is proved from Kakutani’s fixed point theorem. Uniqueness and stability are 
deliberately excluded from this theory, presumable because the former is essentially trivial but messy 
and the latter leads into the vast area of dynamic processes where no orderly system is as yet in sight. 
Chapter 7 treats briefly the optimality of the equilibrium. This is followed in Chapter 8 by an extension 
of the model to a situation of uncertainty which interestingly enough does not involve any probabilistic 
considerations—and hence has no need of the apparatus of expected (cardinal) utility. This final chapter 
gives in effect a condensed summary of the entire theory of this book. 

One can only marvel at the precision, the elegance, and the perfection of this work. Throughout, 
it is an original, penetrating statement of the heart of neoclassical economic theory, the theory of value 
and equilibrium, and it brings out as never before the unity and architectural beauty of this body of 
thought. Its impact will be felt by every mathematical economist. Without question this book will be 
a standard reference. Beyond this it should prove a boon to all mathematically educated persons in 
search of a royal road to economic theory. 

With this latest volume the series of Cowles Foundation monographs—including such classics as 
“Activity Analysis of Allocation and Production” and ‘Statistical Inference in Dynamic Economic 
Models’””—has achieved a new high point. Alfred Cowles, to whom this book is dedicated has indeed 


been honored. 
MARTIN BECKMANN 


(Continued on p. 244) 
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VERTEX EXCITED SURFACE WAVES ON ONE FACE 
OF A RIGHT ANGLED WEDGE* 


BY 
S. N. KARP anp F. C. KARAL, Jr. 
Institute of Mathematical Sciences, New York University 


Abstract. The problem of the propagation of electromagnetic waves by a magnetic 
line dipole source located at the corner of a right angled wedge is considered. It is assumed 
that an impedance or mixed boundary condition is prescribed on one of the wedge 
surfaces and that a homogeneous boundary condition is prescribed on the other. The 
impedance boundary condition is such that surface waves are generated. The amplitude of 
the surface wave generated is determined. A comparison is made between the magnitude 
of the surface wave for this problem and that of a magnetic-line dipole source located 
at the corner of a right angled wedge with the same impedance boundary condition 
prescribed on both* surfaces. The far field amplitude of the radiated electromagnetic 
field is also given as an elementary function of the angle of observation. 

1. Introduction. The problem of the propagation of electromagnetic waves produced 
by a magnetic line dipole source located at the corner of a right angled wedge is con- 
sidered. It is assumed that an impedance or mixed boundary condition is prescribed on 
one of the wedge surfaces and that a homogeneous boundary condition is prescribed 
on the other. The impedance boundary condition is such that surface waves are generated. 
The amplitude of the surface wave is determined. A comparison is made between the 
magnitude of the surface wave for this problem and that of a magnetic-line dipole 
source located at the corner of a right angled wedge with the same impedance boundary 
condition prescribed on both surfaces. The latter problem has already been treated by 
the authors (see Karp and Karal [3]). The comparison of the surface waves for the 
two different configurations is made on the assumption that the sources have the same 
strength. The far field representation of the radiated electromagnetic field is also given. 

The impedance boundary condition prescribed in the present problem is of the form 


1 au 
r 06 


where u is the z component of the magnetic vector, r and @ are the usual polar coordinates, 
a is the wedge angle, and X is a constant characteristic of the surface. The value of \ 


is given by 


—u=0 6=a (1.1) 


A= + tweZ = + twe(R — 1X) (1.2) 


where ¢ is the permittivity of free space, w is the angular frequency and Z, R and X are 
the impedance, resistance and reactance of the surface, respectively. It is well known 
that for homogeneous media of large finite conductivity, R and X are positive in sign, 
small in magnitude and approximately equal. For corrugated or dielectric-coated surfaces, 
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Research Directorate of the Air Force Cambridge Research Center, Air Research and Development 
Command, under Contract No. AF 19(604)1717. 
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R and X are positive in sign and R is much smaller than X. Hence X is positive in either 
case and Eq. (1.2) implies that 
Red > 0. (1.3) 


For illustrations of these boundary conditions see the Introduction of Refs. [3] and [7] 
where the appropriate references are cited. The condition R « X is of course the most 
important special case insofar as surface waves are concerned. 

The problem we treat is not separable because of the mixed boundary conditions. 
This difficulty can be overcome by the introduction of an auxiliary function which is 
a linear combination of the magnetic-field and its cartesian derivatives. The auxiliary 
function is chosen in such a way that it satisfies the wave equation and simple homo- 
geneous boundary conditions. Once the auxiliary function is obtained, the original 
field can be determined by solving a partial differential equation. This idea is due to 
Stoker [6] and Lewy [4] who studied problems in water wave theory. Similar ideas have 
been employed by the authors in the solution of diffraction problems arising in electro- 
magnetic theory. See [1], [2], [3] and [7]. 

Section 2 contains the exact solution of the problem stated in the preceding para- 
graphs. The limiting cases of small and large \/k are studied and simplified expressions 
for the surface waves in both cases are given. A comparison between these results and 
those obtained in the problem of a line source located at the tip of a right angled wedge 
with the same impedance boundary condition prescribed on both surfaces is made. In 
Sec. 3 we discuss the radiated far field and obtain simple expressions for arbitrary values 


of X. 

2. Solution. Consider the right angled wedge y = 0, x > 0 and xz = 0, y < 0, 
and suppose there is a magnetic-line dipole located at the origin x = 0, y = 0. (See 
Fig. 1.) The boundary conditions on the surface are given by 

JU 
Ee 0, y=0, «>0 
oY 
(2.1) 
10u 
lan ~ XU = 9, x= 0, y <0 


/ Line Dipole Source 


pom | 
we eee EP Oe a a a ae al 
Perfect Conductor 
(Homogeneous Boundary 
Condition) 





x 


\ 








Impedance Boundary 
Condition Prescribed on 
Wedge Surface 


FIGURE | 
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where u is the magnetic vector H, and ) is a constant related to the material comprising 
the surface. We wish to solve the time reduced Maxwell’s equations, subject to the 
prescribed conditions, and obtain the amplitude of the surface wave propagated on the 
surface x = 0, y < 0. 

The time dependent form of Maxwell’s equations is 


- H = —wweE, (2.2) 
curl E = iwuH, 
where E and H are the electric and magnetic field intensities, and ¢ and wu are the per- 


mittivity and magnetic permeability of free space. We assume the time dependence 
to be of the form exp (— tw). Because of the geometry, the field produced is independent 


of z and hence the field is completely determined by the value of H, . We have 


H, = H, =E,=0 (2.3) 
and 
. 1 OH, 
E, — a Ou ? 
WWE y (2.4) 
Be, = 4, 
Wwe OX 


The field component H, = uw satisfies the equation 
(V? + ku = —4x d(x — 0) Hy — 0), (2.5) 


where V° is the rectangular Laplacian, k is the propagation constant of free space and 
6 is the Dirac delta function. Therefore, the mathematical problem reduces to that of 
solving the inhomogeneous wave equation (2.5) subject to the boundary conditions 
given by (2.1). In addition to the prescribed boundary conditions, we require that the 
far field be outgoing and that the total magnetic field, excluding the source, be finite 


everywhere. 
Let us make the substitution 


y= OF hu. (2.6) 


Then v satisfies the following equation: 


(V+khy=0, 2 +y #0, 0<0<%, (2.7) 


subject to the simple conditions 
2=0 y=0, 2>0 
, (2.8) 
(2)v=0 z=0, y <0 
(3) Outgoing waves at infinity. 
If the problem for v can be solved, then we can obtain the solution for the problem 


involving wu since (2.6) can be integrated. 
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A particular solution of (2.6) is 


u,(z, y) = —exp (r2) | exp (—Azx)o(E, y) dé. (2.9) 
2 

The range of integration is dictated by the condition that Re \ > 0. It is important to 
point out that the function v does not have to be finite at the origin. Thus we should 
first introduce all radiating solutions, no matter how singular at the origin, that satisfy 
the appropriate boundary conditions for v. Out of this class of solutions we finally select 
only those that yield an everywhere finite value of the total magnetic-field (excluding 
the source) in the original physical problem involving u. By employing the method of 
separation of variables, we find that the most general representation for v which satisfies 
conditions (2.8) is of the form >> a,H‘}.,,)/3 (kr) cos [(2n + 1)/3]0, where n is an integer 
and H‘},,)/3 (kr) is a Hankel function of the first kind of order (2n + 1)/3. Since v 
cannot be too singular, the only admissible functions of this class are 


v(a, y) = coHl{? (kr) cos @ + exHi2A(kr) cos 5 2.10) 


The above expression for v is a wave function, defined, by means of (2.7), in the angular 
sector 0 < @ < 37/2. We note, however, that we can extend the range of definition of 
v to include the angular sector 0 < @ < 2m by means of the reflection principle since v 
vanishes on the negative y-axis. This extension is necessary since the range of integration 
of the particular solution (2.9) extends into the fourth quadrant for negative values of y. 
Using this extended definition of v, the particular solution (2.9) becomes 


/ : oe Ke ( 6 ’ 
u,(x, y) = —exp (dé) [ exp (— nae) cos @ + c,Hj)3(kr) COS 5 dé. 2.11) 

Jz e 
The above result, though well defined, does not satisfy all the conditions of the problem 
since it is not regular across the line y = 0, « < 0. This difficulty can be overcome by 
adding a complementary solution, u, say, of the inhomogeneous equation (2.6) to (2.11). 
We then obtain the representation 


u(z, y) = —exp (Az) | exp (— feo (kr) cos 6 + ¢,H})2(kr) cos 4 dé 
v0 (2.12) 
4. {2 exp [Ax — i(k? +2)'?y], y <0) 
| , 
0 . y > 0) 


where r = ( + y’)'”, 6 = arctan &/y and 0 < @ < 2x. The complementary solution 
u. , Which is indicated in square brackets in (2.12), is obviously not a wave function 
since it is discontinuous across the line y = 0. However, although each solution separately 
is discontinuous together with its y-derivative along the negative z-axis, the solution 
u =u, + u, of (2.12) must be regular there. From this fact we can derive two relations 
between the constants ¢) , c¢,; and c, by calculating the discontinuities of uv and du/dy 
across the negative z-axis and setting the discontinuities equal to zero. These calculations 
completely solve the problem because the constant c, is related to the source strength 
and is known. Hence there are a sufficient number of equations to determine all of the 
unknown constants. We note that the form of the term in brackets in (2.12) is that of 
a surface wave with an amplitude given by ¢, . 

Before proceeding to the calculation of the discontinuities and the value of the 
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constants c, and c, , we first simplify the expression for u and then determine du/dy. 
By using the relation 
© HS (kr) = —kH{” (kr) cos 8 (2.13) 


and integrating by parts, (2.12) becomes 


u(r, y) = -? {15° — \ exp (Az) i: exp (—A£)H4” (kr) as} 


— ¢, exp (Az) ia exp (—At)H})2(kr) cos dé (2.14) 


4 * exp [Az — i(k +2)'7y], yy < i 
0 ‘ y> 0 


At this stage we can determine c) since we recognize the Hankel function Hj" (kr) as 
being the Green’s function for a line source. The solution of (2.5) for a simple point 
source located at the origin is given by the well-known result 


G(r, 0;0, 0) = irH” (kr). (2.15) 
Hence the proper normalization for a point source is obtained provided we let 
Co = —ink. (2.16) 
We next determine du/dy. We have 
Ou, i 2 f ) a) | 
ry=— San eogt Pome a es , on (~~ = 
u(t, y) = 57 ew) i E Hy" (kr) — d exp (Az) j exp (—Aé) ay {Ho '(kr)} dé 


— ¢, exp (Az) / exp (—)é) E {9309 cos 4} dé 
* io + d’)'” exp [Az — i(k? + ’)'*y], oy < | (2.17) 
0 F y > 0 


where we have taken the derivative inside the integral sign. By using the following 


relations 


2 (Ker) = —kHS? (kr) sin 8, (2.18) 

oy 
oe f a 6 PRS . 20 
ay Vi 3(kr) cos a ~_ {1309 sin s} = k exp [i(2x/3)]H2/a(kr) sin 3 (2.19) 
and integrating by parts we obtain 
u(r, y) = 3 { — RH (kr) sin 6 + kd exp (Az) [ exp (—A£)H;” (kr) sin oat} 

(1)  - - (1) iw 
+ c,H,/3(kr) sin - CA exp (Az) [ exp (—A£)H1/3(kr) sins dé 


(2.20) 
° | 2 
— ek exp [i(2x/3)] exp (Az) / exp (—E)HS23(kr) sin = dt 
| eat(h? + 0°)? exp [Ar — i(k? +°)'7y], oy < 4 
| 0 ' y>oO 


= 
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It was noted earlier that the constants c, and c, can be determined by requir- 
ing that u and its derivatives be continuous throughout the region of physical space 
(0 < @ < 3n/2). This is equivalent to calculating the discontinuities of u and du/dy 
across the negative z-axis and setting the discontinuities equal to zero. Hence for u(z, y) 
to be a wave function, we require that the following jump conditions be satisfied: 


| (I)[u] = {u(z,0*) —u(@z,0)} =0, 2<0 (2.21) 
(II) fu,] = {u,(a, 0°) — u(x, 0°-)} = 0, x<0. 
The first jump condition yields immediately 
Cc. = —$e,1,,3(k, d), (2.22) 
where we have used the notation 
I(k,d) = [exp (—d)HSk | € |) a. 2.23) 


We note that 

eras t exp [—7r(2/2) 
1k, ) = <i exp L=ile/2) 

sin vr k’(\? + k?)'? (2.24) 
- {exp [—iv(r/2) (IX + ]'?— d\)’— exp [iv(r/2)]((A? + k7]'? + A)’, 

where —1 > Rev < + 1. (See Magnus and Oberhettinger [5].) The second jump con- 
dition requires more careful treatment. If we substitute (2.20) into jump condition 
(II) and use (2.16) we obtain 


= Lim 277k) exp (Az) . xp (—r\S/) ASP (ke é fe — £ 
[u,] im 27ikd exp (Az) | exp (—AE)H, (K[E + y'] €+ yy” dé 


y+0 dz 


oe T seated 
+ ¢, exp (Az) sin 3 | exp (—A£)Hy/(k | & |) dé 





2.25) 
,. a 2” 
+ ¢,k exp [i(2x/3)] exp (Az) sin 7 | exp (—A£)Hy)3(k | & |) dé 
+ coi(k? + d*)'” exp (Az) = 0. 
Let us consider the expression 
= Li [ x (—né Yk P oie > Y 2 > ‘2 96 
J im exp (—Af)A, (RIE + J”) f+ 7)" dé. (2.26) 


yotO0 dz 


The integrand in (2.26) tends to zero with y except in the neighborhood of § = 0. Hence 
the only possible contribution occurs when £ lies in an interval — e to e, where « is arbi- 
trarily small. Since we assume that y is already small, we can use the asymptotic repre- 
sentations for the exponential and Hankel functions appearing in the integrand. Hence 


€ 


My; , y 
=-— 3 OE (¢ 7 
J 7 Lim [ Pay dé. 2.27) 


If we make the substitution = yl, the above integral becomes 


y7+0 


af t 
= —_— ee 9 92 
J 21.747 i (2.28) 
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Substituting this result into (2.25) and using the notation defined in (2.23), the second 
jump condition becomes 


9 
+ 4rd + ¢, sin 5 I, ;.(k, ®) 
(2.29) 


+ ¢,k exp [i(22/3)] sin = Toja(k, +) + coi(h’ + *)'? = 0. 


3y substituting (2.22) into (2.29) we find c, . We obtain 











c = aa. ES TT ae ee 2 (2 30) 
a 3 (A — 1k + NY) Lis — & exp [i(24/3) o/s ve 
and 
9 
ia (2.31) 


i + 2 31 9 9172 -77,2 23172 . ) a 
A— 3k +r) )lis — k exp [t(2r/3) o/s 


The value of c. is the amplitude of the surface wave propagated along the wedge surface. 
We next determine the values of c, and c, for small and large values of A/k. When 
\/k is small it is easily shown that 


fiat : (1 — 13°”) = 23° exp [—i(n/6)] (2.32) 
1 *21/2 2 . 
Toj3 FZ i (il — 73°) = j, oP [—i(r/3)] (2.33) 
and the approximate expressions for c, and c, become 
C, &Y an exp [—7(x/3)] (2.34) 
Co & 2ri(>), (2.35) 
When i,k is large we find that 
— es oe 
I, so — $s 3i72 (?) ? 3 (2.36) 
=F: we " 
Ii, — 1 3172 (2) ys (2.37) 
and the approximate expressions for c, and c, become 
2/3 
¢, &e wk exp (—ite/)(*) (2.38) 
C, &Y 2 3'r exp [i(/3)]. (2.39) 


It is of interest to compare the amplitude of the surface wave generated in this 
problem with that generated by a line source located at the tip of a right-angled wedge 
with the same impedance boundary condition prescribed on both surfaces. The latter 
problem has previously been solved [3] and the amplitude of the surface wave propagated 
on the wedge surface x = 0, y < 0 is given by 
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1/2 


C3 = 4: arcs are 7; : per r. a ) 
3 {[a(k’ + d*)'? — d3°* ]I2,3 + 3'k exp [i(r/3) ]1,3} (2.40) 
. a ae 
(i(k? + r7)'* — dJ 


Fas Al o/s a k exp (t(r/3) Wiss 


When }/k is small we find that 


Cs & 2ri(*) (2.41) 


L 


and when }/k is large 


Ar (3)" a 
eS = 2 exp [2(57/12)]. 2.42 
(1 — a)(1 — 73°”) "\2 P [s(5x/12)] ) 
Hence, when the value of \/k is small, the ratio of the amplitude of the surface wave 


on a wedge with an impedance boundary condition on one surface to that on both 
surfaces is given by unity. When the value of \/k is large, the ratio is given by 2°’~. In 
the limit \/k — © therefore, the total energy that goes into surface waves in our present 
configuration is the same as the energy so utilized when both faces of the wedge can 
support surface waves. 

3. Determination of the radiated far field. We have already obtained the exact 
solution for a line source located at the tip of a right-angled wedge with an impedance 
type boundary condition prescribed on one surface. The exact solution as given by (2.14), 
however, is complicated and it is desirable to find a simple expression for the radiated 
far field. From (2.6) and (2.10) we have 


yu , 
v= — — hu 3.1) 
Ox 
and 
(1)/ | 1) \ 6 9 
v = ¢H;'’(kr) cos 6 + ¢,H;/3(kr) cos 3" 3.2) 
When kr is large (3.2) becomes 
9 1/2 
iiss (2) exp {i[kr — (3m/4)]} cos 8 
TK 
3.3) 
2 ; 6 
+ ¢,|— exp {i[kr — (52/12)]} cos=- 
akr o 
We expect the far field of the z-component of the magnetic field vector to have a similar 


form, except that it contains surface wave terms which are appreciable near the wedge 
surface x = 0, y < 0. Hence for a fixed 0, as r — ~ we have 


u= E exp (tkr) + o(,) | 6 # 3n/2, 3.4) 
where /(@) is unknown. Now 


0 0 ts 0 —— 
_* cos 6 — 7, sin 0-5 (3.5) 
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If we substitute (3.3), (3.4) and (3.5) into (3.1), and keep terms of order r~’”’, we obtain 





(3.6) 


1/2 Co exp [—i(32r/4)] cos 6 + c, exp [—7(52/12)] cos 
) ik cos 6 — X 


Hence the desired radiated field is given by 





: , 6 
9 \1/2 +o exp [—1(37/4)] cos 6 + ¢, exp [—7(5r/12)] COs 3 1 
ulr. 6) = = e* - : +0 a7 2 (3.7) 
kr ik cos 6 — X r 
where ¢y . ¢; and ¢c, are given by (2.16), (2.30) and (2.31), respectively. 
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Applications of the theory of matrices. By F. R. Gantmacher. Interscience Publishers, 
New York, London, 1959. ix + 317. $9.00. 


This is a translation of the second part of the author’s A Theory of Matrices which was published 
in 1954 in Russian. As such it constitutes a unique mathematical monograph. The applications alluded 
to in the title lie in other fields of mathematics. The reader should possess considerable maturity in 
these fields as well as in matrix theory in order to appreciate the elegance of the author’s developments. 

There are five chapters. The first two, treating the normal forms of complex matrices and singular 
pencils of matrices, are quite formal being little more than a sequence of lemmas and theorems. The 
third chapter treating matrices with non negative elements begins formally but becomes more readable 
when it turns to applications to Markov chains with a finite number of states and to oscillating matrices. 
The two final chapters contain the major applications. Chapter 4 is a study of systems of linear differ- 
ential equations and Chapter 5 is devoted to the Routh-Hurwitz problem. 

As an indication of the level, the chapter on differential equations includes Lyapunov’s stability 
theory, the theory of singular points and works up to conclude that certain results of V. Volterra and 
G. D. Birkhoff concerning isolated singular points were erroneous due to insufficient examination of the 
degenerate cases involving nonlinear divisors. 

The author’s developments are elegant and worthy of study by experts in the fields mentioned. 
In one case they appeared overly elaborate to this reviewer. In Chapter 3 the properties of oscillating 
matrices are developed to show that certain vibratory systems (e.g., a stretched string carrying n masses) 
possess the Sturm properties. Apparently the author was unaware that in 1884 E. J. Routh had demon- 
strated these same results in a much more direct manner by extending Sturm’s argument to finite 


difference equations. 
STEPHEN H. CRANDALL 


Matrix calculus. By. E. Bodewig, Second Edition. North Holland Publishing Co., 
Amsterdam, and Interscience Publishers, Inc., New York, 1959. xi + 452 pp. $9.50. 


Recent progress in digital computer techniques for the determination of the eigenvalues and -vectors 
of large matrices is reflected in this new edition by the inclusion of Lanczos’s pq-algorithm, Rutishauser’s 
LR-algorithm and Wilkinson’s method. Several hitherto unpublished results by the author are also 


described. The work retains its character as an invaluable fund of information on matrix methods and 


their relative merits in different context. 
WALTER FREIBERGER 


Mathematical methods and theory in games, programming, and economics. By Samuel 
Karlin. Volumes I and II. Addison-Wesley Publishing Co., Inc., Reading, Mass., 
and London, 1959. x + 433 pp. (Vol. I), xi + 386 pp. (Vol. IT). $12.50 each. 


The last twenty years have seen the remarkable growth of a new branch of applied mathematics— 
the use of mathematics in decision making within both business and government. These two books 
present a systematic and elegant account of many of the mathematical methods which have been de- 
veloped for this subject, written by one of its leading practitioners. 

The first four chapters of Volume I are fairly standard material about matrix games. Then come 
linear programming—theory and computational methods—and a chapter on nonlinear programming. 
The various famous special types of problems (optimal assignment, transportation, network flow, etc.) 
are taken up in detail. Fenchel’s theory of conjugate convex functions is developed from a new geo- 
metric viewpoint, and is applied to get a duality theorem for nonlinear programming. The rest of Volume 
I deals with theories of production and consumption, with welfare economics, and with the dynamics 


(Continued on p. 270) 
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A GENERAL THEOREM CONCERNING THE STABILITY 
OF A PARTICULAR NON-NEWTONIAN FLUID* 


BY 
SAMUEL M. GENENSKY 
Mathematics Division, The RAND Corporation 


Summary. It is the intention of the present paper to prove a theorem concerning 
the stability of a particular non-Newtonian fluid suggested to the author by Professor 
R. 8S. Rivlin of Brown University. The method used in proving this theorem is similar 
to that employed by H. Schlichting in his proof of a similar theorem for an inviscid 
fluid which was originally established by Lord Rayleigh. The acceleration gradients 
introduced by the non-Newtonian fluid model into the constitutive equations are found 
to alter the stability criterion set forth by Rayleigh for an inviscid fluid. 

1. Introduction. As early as 1880 Rayleigh [1] proved that for an inviscid fluid 
the existence of a point of inflection in the velocity profile of a steady one-dimensional 
basic flow is a necessary condition for the growth of a superimposed two-dimensional 
disturbance. It is the intention of the present paper to prove a similar theorem for a 
particular non-Newtonian fluid suggested to the author by Professor R. S. Rivlin.* The 
method used in proving this theorem is similar to that employed by Schlichting [2] in 
his proof of the Rayleigh theorem. 

2. The constitutive equations and equations of motion. Let X,(j = 1, 2, 3) be 
the coordinates referred to a rectangular Cartesian coordinate system 2; , of a generic 
particle of a continuous medium in the undeformed state at time f, . Let x; be the co- 
ordinates of the same particle in the deformed state at time ¢. It then follows that the 
components of velocity v; and the components of acceleration a; of the particle are 


given by 
Ox; 0°x,; 
vo, = — and 4a, = 7, (2.1) 
ot ot 
where xz; are regarded as single-value continuous functions of X,(k = 1, 2, 3) and @, 


having as many continuous derivatives as the analysis requires. As is well known, if 
v; are considered to be functions of x, and ¢, then 


Ov ov 
= — ) —. (2.2 
a; at " OX, ) 
Now the equations of motion are 
Ov; ov \ Ot;1 / 
oe f, = —it l= 1,2,3 (2.3 
(% + taz,) + Pb ax,’ ( ’ ), \ ) 


where p is the mass of the medium per unit volume and f; are the components in the 
coordinate directions of the body force per unit mass, also measured in the deformed 
state. The components of stress ¢;, resulting from the deformation are defined as follows; 
t;, , tj and ¢;, are the components of the force per unit area in the positive direction of 
the xz, , 2, and x, axes respectively, measured in the deformed state, exerted across an 





‘Received October 28, 1958; revised manuscript received August 12, 1959. 
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element of area at (x, , X2 , X;) normal to the 2; axis, by the material on the positive side 
of the element upon the material on the negative side of the element. 

Rivlin [3] showed that if ¢;, at the point x, and at time ¢ are assumed to be poly- 
nomials in the velocity gradients dv,,/dx, (m,n = 1, 2, 3) and the acceleration gradients 
da,,/dx, and if, in addition, the medium is assumed to be isotropic at time /) , then the 
stress matrix T = || #,, || is expressible in the form 
T= gI+ giA, + PoAs + ¢3Aj + gsAz + gs5(A, A, + A.A,) (2 1) 

2 2) i 2 2 2,2 2,2 : 
+ vo(Ai A, + A.A)) + ¢7(A, A; + A;A)) + gs(A, Az + A,A)), 

where I is the unit matrix, A, and A, are symmetric kinematic matrices defined by 
ll ¢ ‘ ; | 
| 0a; Oa, > Wm Wm a 
A, = || = + mm Ae +S +3 (2.5) 

1] 02x, On, | || Ox, Ox; Ox; O02, | 
and ¢, (¢ = 0, 1,2, --- , 8) are polynomials in ¢rA, , frA, , trAj , trA?, trAj , trA; , trA,A,, 
trA{A, , frA,A? and trA{A; . In a later paper Rivlin [4] pointed out that for an incom- 
pressible material, the stress corresponding to a specific state of flow is indeterminate 
to the extent of an arbitrary hydrostatic pressure p. Since in the present paper we shall 
confine our analysis to an incompressible fluid, we may replace gp in equation (2.4) by 
— p. We shall also restrict our investigation to a fluid for which ¢g, and ¢, are constants 
and ¢, (q = 3, 4, 5, --+ , 8) are identically zero. 

Equation (2.4) now takes the form 


T = —pI + ¢,A,; + gA: , (2.6) 
or alternatively 
Ov; Ov; 0a; Oa, OV», Om = 
i, = —pé + on( 2 4 2) + (28 4 § egies ; (2.7) 
i sili *N\ Ox, Ox, *\ 02; Ox; OX; OX, 
Introducing Eqs. (2.2) into Eqs. (2.7), we obtain 
‘a Ov; , Ov, Ov; , Wm Ov, 
n= =? in +o 21 +e += 
OX Ox; Ot OX, OL; OL» (2.8) 
av, 0°v, OV, OV Ov, <> Um Wm 
+ 0m i a pt te + 2). 
O02; 02, Ot Ox; OL; OX m OX; CLia OX; OX, 


Since we have assumed that the fluid is incompressible, the continuity equation is 


4 
=" = 0. (2.9) 
ie 

Introducing the constitutive equations (2.8) in the equations of motion (2.3), employ- 


ing the incompressibility condition (2.9) and neglecting body forces, we obtain 


(2 = 21) Op . Ov ( O’v, 
PL TRS Se tee ee oe 
at Ox, Ox * Ox, Ox; F\ Ot Ax, OX, 
— WW; 8m g Wm Ov; 
ae er a ee 
IX m OL, OX OX, OX, OX 
sii lade ici, ie ei (2.10) 
he re) y ; Vm 0 *v, 
Um > ye a? ‘ 
OX, OX, OXm OX, OX; OX», 
492 Wm 9m Wn ml ) 
Ox, OX; OX, Ox; OX, OX, 
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3. Development of the stability equation. Consider next a two-dimensional steady 


laminar flow with velocity components 


W, = W,(z.) and W.= W; = 0. (3.1) 


Examples of such flows are the flow between a pair of parallel plates sufficiently removed 
from the intake section, and the flow in the boundary layer along a flat plate excluding 
the region of its leading edge. The pressure P necessary to maintain such flows is in 


general a function of both x, and z, . 
We shall superimpose upon the laminar flow a small two dimensional disturbance 


with velocity components and associated pressure given by 


Uy u(t, , 2,8), Us = U(x, , 22, 8); uz; = 0 and p* = p*(z, , 22,4). (8.2) 
Thus the velocity components and pressure of the composite flow are 
n=Wit+u, Ve = Us, v,=0 and p=P+>7*. (3.3) 


Further, we shall require that the composite flow satisfy the same boundary con- 
ditions as the steady laminar flow. Thus the disturbance satisfies the boundary conditions 


u.=u=O0 at 2 =0 and z=L (3.4) 
in the case of flow between paraliel plates separated by a distance L, and 
u.=u=0 at z%=0 and wz 0 (3.5) 


in the case of boundary layer flow along a flat plate. 
Next, introducing the velocity components and pressure given by (3.3) into Eqs. 
2.9) and (2.10), and assuming that the velocity components wu, are sufficiently small 


so that terms of the second degree in u; and derivatives of u; may be neglected in com- 


parison with terms of the first degree, we obtain 
ou OU Qf 
L +. SS ee @, (3.6) 


OX, OX 


a on AP 
(2 + W, Se + Win) = Oe 4 tw + WO 


Ot OX, OX, 





‘ e( 7? % 2 Oi 4 owe % 4 aw? rs 
lene (3.7) 
+ WV + 30! oe ‘+ WV? 2 Si 
vy 
+ Wf’. + 2WI 4) 
Ox; 
ind 
OUp . - OU» \ _ _ oP _ op* ( 2 OUs 
o( Ot - n ; a, OX OX + av “+ sls Vv ot 
4 4! as — + sw +W, ‘i (3.8) 
v 





+ 4WiWy’ + aa ot ns 2W! £ 2 + 2Wiv* ‘u), 


2 
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where VY’ 
tO 2%. 
Further, 


laminar flow 
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(9°/dx3) + (d°/dx;), and primes denote ordinary differentiation with respect 


introducing Eq. (3.6) into Eqs. (3.7) and (3.8) and assuming the basic 
satisfies the equations of motion, we obtain 





Ou, es Op* 2 OU, 
o( 2 ra +e us) = oe +aVu + | Vo = 
a Or; ‘ OX, ot = 
(3.9) 
— OU; ppp Clee r 72 ON, , —- —— 
+ 3Wi —— + Ww! — H+ WV t+ OW ii, + BW 3 
OX, 0X2 OX; OX, OX} 
and 
OU> _ ou op* 2 [ 9 Othe ; Oru 
om + W, — 9 ~~ g,V Us T P2\ PY, 2 ar. 3d = 
dl OX 7>1/ OX2 dl OX, OX2 . ‘ 
: shies (3.10) 
Le a | a ug Gee 4 Ou, 
+ 3Wi’ — + W, —e oe AWS’ — + 4W 1 — ) 
OX, OX OX2 OX> 
Differentiating Eq. (3.9) with respect to x, and Eq. (3.10) with respect to x, , elimi- 
nating 0°p*/dx,dx, from the resulting equations and again employing Eq. (3.6), we 


arrive at 


0 a re] a O 7 . @ } Ou OUs . . 

2 i Hy + 2 )( one 2) = ( Ww — Wu, (3-11 
(2 Ge . 7V ot y ; O21 OX Ox, Y U2 ) 
where v = ¢,/p andy = ¢/p. 


Further, we shall assume that the velocity components of the disturbance are of 
the form 
, 4) = u¥(z.) exp [7A(x, — Cd], (7 = 1, 2) (3.12) 


U4, , Xe 


where | u*u* | is the amplitude and A the wave number of the disturbance, and C is a 
complex number. The real part of C, C, , is the phase velocity of the disturbance and 
the imaginary part, C; , is the amplification factor. If C; > 0, the disturbance tends to 
grow; if C; < 0, the disturbance decays; and if C; = 0, the disturbance is neutral. 
Introducing the velocity components (3.12) into Eqs. (3.6) and (3.11), we obtain 


me ae SP (3.13) 
and 
tA(W, — C)(u*’ — iAu*) 
— + iAYW, — C)\(ut?” — A*ut! — iAut” + iA’ut) (8.14) 
= (yW,"” Wi us. 
Replacing u¥ in Eq. (3.14) by tux’/A from Eq. (3.13), we obtain 
tA(W, — C)uz’”” — A°u§) — [» + tAy(W, — C))uk""” — 2A°uz” + Ave) (3.15) 
LA(Wh! — yWh"” ui 
Writing Eq. (3.15) in dimensionless form by letting - nfl vy, W./W, . 


where L and W, are a characteristic length and 


a 


AL and c C/W 
velocity respectively of the steady laminar flow, we arrive at 


We W 
a characteristic 


, & 
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ia(V, — c)(wh’ — a’w.) — E + = (V,; - 0) [ows — 2a*wi’ + a‘we) 
(3.16) 


~ ia( V1 -3 v5", 


where R = W,L/» is the Reynolds’ number of the laminar flow, S = L’/y is a non- 
dimensional parameter arising from the presence of non-Newtonian terms in the con- 
stitutive equations of the fluid, and primes denote ordinary differentiation with respect 


tO No. 
We see that if S — , Eq. (3.16) becomes the familiar Orr-Sommerfeld stability 


equation. 
Finally, in terms of w, the boundary conditions (3.4) and (3.5) become 
w.=w,=0 at » =0 and » = 1, (3.17) 
and 
w2=w,=0 at m=0 and m->@ (3.18) 


respectively. 
4. The general theorem. Under the assumption that S is finite and RF is infinite, 


Eq. (3.16) takes the form 


, f , 2o4/ yer 1 , 
(V, —o(wi’ — aw) - ‘ (V; —c)(wi’” — 2a’w}’ + a‘w.) = (Vi! - s Vi"'w,. (4.1) 


We shall now prove the following theorem. 


The existence of a point in the flow field for which V’) — (1/S)V{’"’ is equal to zero, 
is a necessary condition for the amplification of a disturbance. 


M(w-) < ww,’ — (1 + 2a) ot + (.° + a + Aer Wa (4.2) 


‘ (vy ae ‘ vy) 


oo _ I atrth 2a” <tr 2 “i . 
M(w2) = % Ws (1 + =¢ )as + (. + ~G)i. + eae De ; (4.3) 


Regarding w, as a complex variable, we define 





Il 


and 





where a bar denotes the complex conjugate of the corresponding unbarred quantity. It 
is easily seen from Eq. (4.1) that both M(w.) and M(w.) are equal to zero. 

To prove the theorem we assume that V’/ — (1/S)V{’’’ ¥ 0 throughout the flow 
field. Then, multiplying M(w,) by @, and M(w,) by w, and subtracting the resulting 


expressions, we obtain 


=, ee in 
w,M(w.) — w.M(w.) = re (w.w3’’’ — we0s’"’) 


2 2 
- (1 a 2a) a. — Wztb2’) (4.4) 


ree 1 rrr) ] in 1 ) ». |? 
+(¥ g "Ni —2 7 V, — a) | T- 














250 SAMUEL M. GENENSKY [Vol. XVIII, No. 3 


Integrating Eq. (4.4) with respect to 7. between the limits 7, = 0 and 7, = 1, we 
have 
1 5 1 {1 
[ f.M(w.) — w.M(w.)] dn. = < [(dow3’’ — wes!’) — (ww! — whd’)] | 
A 0 


( rer 1 vir") 
al 1 S 1 


i 2a” - , -/ [* ° 2 
-—i oe (wow, — ws] , + 2ic; | ee APS ee | we | dye . 


70 


(4.5) 


A 


Because of the boundary conditions (3.17), the first two terms on the right hand side of 
Eq. (4.5) vanish, and further, since both M(w.) and M(w.) are equal to zero, the left 
hand side of the equation is equal to zero. It then follows that 


(vy ae ‘ yi") 
= : | 


1 
C; [ ss , ee {2 ea | We ‘ dno 
|Vi-—e 


v0 


must vanish. But this is impossible, since c; > 0, for a disturbance which tends to grow, 
and hence [1/| V,; — c |*] > 0. Further | w, | is positive and by assumption V’/ — 
(1/S)V{’’’ # 0 everywhere in the flow field. Thus we conclude that for a disturbance 
which tends to grow, there exists an 72,0 < n. < 1, for which V’// — (1/S)V{’”’ is equal 
to zero. 

In the case of laminar flow by a flat plate the proof of the theorem differ only in that 
the limits of integration become 7, = 0 and 7, = ©, and the boundary conditions 
(3.18) are used in place of those given by (3.17). 


{EFERENCES 

1. Lord Rayleigh, On the stability or instability of certain fluid motions, I, Proc. Lond. Math. Soc. 9, 
57-70 (1880); see also Scientific Papers 1, 474-487 (1899) 

2. H. Schlichting, Boundary layer theory, translated by J. Kestin, McGraw-Hill, Inc., New York, 1955, 
pp. 314-323 

3. R.S. Rivlin, Further remarks on the stress—deformation relations for isotropic materials, J. Ratl. Mech. 
and Anal. 4, 681-702 (1955) 

4. R.S. Rivlin, Solutions of some problems in the exact theory of visco-elasticity, J. Ratl. Mech. and Anal. 
5, 179-188 (1956) 


251 


TRANSIENT MOTION OF A LINE LOAD ON THE SURFACE 
OF AN ELASTIC HALF-SPACE* 


BY 
DANG DINH ANG 
Guggenheim Aeronautical Laboratory, California Institute of Technology 


Abstract. The present paper studies the wave patterns generated in an elastic 
half-space by a line load moving on its surface with a velocity varying as a step function 
of time. The solution given in closed form is obtained by means of Fourier integral 
equations techniques following a Laplace transformation with respect to the time 
variable. The inversion of the Laplace transforms is based on a trick due to Cagniard 
and De Hoop. 

Introduction. For the past few years, people have been interested in the effects 
that moving blast waves on the surface of the earth exert on its interior and in particular 
on underground structures. The problem has been considered by Sneddon [1] and 
independently by Cole and Huth [2], who treated the steady motion of a line load on 
the surface of an elastic half-space. Although it is clear that, far enough from the starting 
position of the load, the motion can be considered as essentially steady, there is one 
important feature that is not exhibited by a steady state solution, i.e., the “resonance” 
effect present in a motion at a velocity approaching the Rayleigh wave velocity. The 
present paper will take account of the latter effect by considering the transient motion 
of a line load at a velocity varying as a step function of the time. 

The initial and boundary value problem. Let a Cartesian system of coordinate. 
x, y, z be defined such that the elastic half-space is represented by y > 0, that the axis of 
the line load is in the z-direction and that its initial position coincides with the z-axis. 
The strength of the load is assumed constant so that the problem is a two-dimensional 
one in x, y. Furthermore, since the medium extends indefinitely in the z-direction, the 
problem is one of plane strain. Then, if « and v denote the components of the displace- 
ment in the z- and y-directions respectively, it can be shown [3] that the equations of 
elastic motion are satisfied if 


_ dO® ow ~~ om OV (1) 
Ox oy’ oy ox’ 


y= 


where ® and W are solutions of 


2, _ 29 ® o- 398 





In (2), v; and v3 are defined as 
vi = p/(X+ 2), 2 = p/n, (3) 


where p is the density of the material, and \, » its Lamé constants. 


*Received March 10, 1959; revised manuscript received August 28, 1959. This investigation was 
supported in part by the United States Navy Bureau of Ordnance through the U. S. Naval Ordnance 
Test Station, China Lake, Calif. 
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The initial and boundary conditions are: 


a = — d(x ~ co ty = 0,2 > @ (4) 
=0 y=0,2 <0 
Tzy = 0 y =0, forall z, (5) 


where 7,, and 7,, are respectively the normal and shear stresses at the surfaces of con- 
stant y, 5 is a Dirac delta-function, and c;" is the velocity of the moving load assumed to 
be a constant. We thus consider a normal load only, but it will be seen that the general 
case could be handled just as easily, and that the characteristic features of the wave 
patterns of direct interest here are not changed if we add to the load a tangential com- 
ponent. To the conditions (4)—(5), we must, of course, add the condition that the waves 
be outgoing waves. To complete the formulation of the problem, we require the wave 
functions ® and W to be continuously differentiable at least twice everywhere except 
at the the wave fronts, and the displacement to be a continuous function of the space 
coordinates everywhere inside the solid. 

The method of solution consists in applying a Laplace transformation to the wave 
equations (2), and then expressing the solutions of the resulting equations as super- 
positions of “‘plane waves” with the amplitude spectra as the unknowns. Thus, under a 
Laplace transformation defined in the usual manner, the Eqs. (2) become: 


V’*é* = vip’s*, V'u* = vip v*, (6) 


where p= is the parameter of the transformation and the superscript* indicates a Laplace 
transform. Let the solutions of (6) be of the forms: 


$* = [ P(s) exp {—p(s’ + v3)'’y + ipsx} ds, (7) 


y= / Q(s) exp {—p(s? + 02)"”y + psx} ds, (8) 


where P(s) and Q(s) are unknown functions and the paths of integration are along the 
real axis. The functions (s* + v7,,)'”” are given non-negative real parts along the path 
of integration; this can be achieved by defining their branch cuts as follows: Re s = 0, 
v,,.2 < |Ims| < © respectively. Finally, the parameter p in (7)-(8) is taken to be a 
real and positive quantity of magnitude sufficiently large to insure the convergence of 
the integrals. 

The problem of course is to determine the functions P(s) and Q(s) from the boundary 
conditions. Now, the conditions of (4) and (5) can be expressed in terms of ® and V 
using (1) and the well-known strain-displacement and stress-strain relations. Assuming 
that (7) and (8) can be differentiated under the integral signs, we accordingly obtain 
for the Laplace transforms of (4)—(5) 


[ {(s? + v3/2)P(s) — is(s’ + v3)'Q(s)} exp (ipsx) ds 
ws (9) 


(cov2/2pp’) exp (—peot) for xz > ¥ 
= 0 z<0 
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[ ” {Qis(s? + 02)*P(s) + (2s? + 02)Q(6)} exp (ipsz) ds = 0 (10) 


the latter equation holding for all x. Equation (10) is evidently satisfied if the integrand 
is identically zero, i.e., if 


P(s) = (8° + 02/2)R(8), — Q() = -—is(F' + v})'? RO), (11) 
where R(s) is still largely arbitrary. Substituting for P(s) and Q(s) into (9) gives: 


F F(s)R(s) exp (ipsx) ds = (A,/p’) exp(—pooxr), 2«x>0 (12) 
= 0 zx<0 
where 
Ay = Cot2/2p, = F(8) = (8° + 03/2)? — 8°(8 + 03)'7 + 02)”. (13) 
The solution of (12) is immediate 
R(s) = (Ao/2mip*){(s — te)F(s)}~*. (14) 


From (11) and (14), the expressions for P(s) and Q(s) can be deduced and the problem 
is formally solved. It is not difficult to show that this solution converges and is the 
actual solution of the boundary value problem. 

The stress wave patterns. It remains to find the inverse transform of the solutions. 
The Laplace transforms of the stresses are given in terms of integrals of the form 


[* = f K(s) exp {—p(s’ + v°)'?y + ipsz} ds, (15) 


where v = v,,2 and K(s) is some algebraic function of s only. The present inversion of 
(15) back to the é-plane is based on a trick originally due to Cagniard [4] and modified 
by De Hoop [5]. The trick consists in so deforming the path of integration that J* is 


reduced to the form 
J? = [ I(t; x, ye dt, (16) 
0 


where / is real and J(t; x, y) is some suitable function. The inverse of J* can thus be 
obtained by inspection. In the course of the deformation of the path of integration, it is 
essential to know all the singularities of K(s). These can be recognized by inspection, 
with the exception of the poles s = + is zg , Sg > v2, + is z being the zeros* of F(s) 
defined in (13). Once the singularities of K(s) are known, it only remains to follow steps 
quite similar to those taken in [6] where the inverse transforms of integrals of the form 
(15) have been found. We refer the reader to [6] for these details and are content to 
state here the final expressions for the stresses corresponding to a velocity of the load 
smaller** than both velocities of sound in the solid 


*Thus sg™ is the well-known Rayleigh wave velocity. 
**The case of a load velocity greater than either or both velocities of sound could be handled just 
as easily but is omitted for simplicity. 














254 DANG DINH ANG [Vol. XVI 
( : 
Co } ¢gin 6 M,(s 
ty, = ——HAti- vr) 2 = THis Rez a 4 
a 7 \(t — vy") (co + is,)(Se + is;) 
M s 
= 1(8,) Cos os — 
(co * is,)(Sp re 18;) 
} S 6 M, M.(s, COs 0 
Ar = 22 H(t — vg) oe ; Re em Hag ee 
mr \(t — vor’) (co + i82)(Sp + is») (Co + %82)(8e + 2 
Cc tsin 6 Mw 
a 7 aa s 07-Im 
: nOer 5 Py? 3 + cos \ n i 23 io. 
eek | ¢sin 0 N,(s,) 
T = (t{—1 ) THD Re- > eee ee 
oa ur H " ) t — vir)” (Co + 18;)(Se + 181, 
« Rate N i(s,) cos 6 } ee 
(Co + 18; Sr + is,) 
tsin 6 N .(s,) 
= — t — v,7r) ee } 7) er ae 
Tey ime H 17 1 > — vr)? Re a aE 
N,(8,) cos 
— Im7 a yf 
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; tsin 6 Nw.) 
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TT J N\+Gn —ee ve J (Co + tw. Se + 1W . 
where 
r=(¢° +7)", @é=tan'(y/z), O<6<z, 
S20 = (/r° — v;.2°)'”” sin 6 + i(t/r) cos 6, 
M,(s) = (8° + 0;/2)’(8e + is)/F(s), 
F(s) being defined in (13) 
Ms 2 + v1)'7(8? + 03)'"(8e + is)/F(s), 
+ vi .)'? = {tsin 6+ 1(f — v}.2r°)'” cos 6} /r, 
we ((A2_ 4 Rp? jie A. .)1/2 2 2 ye \1 
2 + 03 2 2 
| (t’/r’)(sin’ 6 — cos’ 6) + v3., — vj.2s8in’ 8, 
B, (t°/r? — v;.2)'"*(t/r) sin 28, 


the + signs in (25) corresponding to cos @ > 0 and cos 6 < 0 respectively. 


II 


{ 





No. 3 


i - 


(17) 


(18) 


19) 
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(f) =1lts <t<or 


ts = vr | cos 6| + ( — v})'“rsin 6 (26) 


Il 


= 0 otherwise 


| 9(8) =-10<0<0;s or r—O03;<0< 2 
6; = cos’ (v,/ve) (27) 


= (0 otherwise 
we = ifF} — /r*)'” sin 6 + (t/r) cos 6}, (28) 


w. corresponding to 0 < 6 < @s and zr — Os < 6 < = respectively. 


ws + vi)? = +1(|w.  — 0)", (29) 
(wl + 02)'? = &% — | we |"), (30) 
N,(s) = (8° + v1 — 02/2)(8° + 02/2)(6e + i8)/F(s), (31) 
N.(s) = is(s’ + ¥})'7(s° + 03/2)(8e + 18)/F(s). (32) 


In the above relations the square roots of real quantities are taken as positive. 

From (17)-(19), it is clear that the wave patterns consist of two cylindrical waves 
radiating from the origin r = 0 with velocities v;' and v;* respectively, two head waves 
whose fronts are tangential to the shear cylindrical wave—the head waves are represented 
e.g. in Ar in (18)—the two surface waves propagating in the positive and negative 
z-directions respectively. The complete picture is shown in the figure below. 





| + ——v' it 
[ R <-Sot-n oR 
4 + x 


il Sa t——| mu 











y 


Stress wave patterns. I. Shear wave, II. Dilatation wave, III. Head waves, R. Rayleigh wave. 


It can be verified that the stresses are of the order (vr? — #*)~’”, for ur > t,v = v4.2. 
This singularity is, of course, integrable and does not affect the continuity of the dis- 
placements. The resonance effect anticipated in the Introduction results from the 
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presence of the factors (cy + is)~* (sg + is)~’, s = 81,2 , in the expressions (17)-(19) 
since — isp is not a zero of the functions M(s), N(s) defined above. 
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THE STOKES FLOW ABOUT A SPINDLE* 


BY 
W. H. PELL anp L. E. PAYNE** 
(National Bureau of Standards, Washington, D. C.) 


I. Statement of the problem for an axially symmetric body. The Stokes flow of 
a viscous, incompressible fluid about a body is defined by the assumption that inertial 
effects are negligible in comparison with those of viscosity, or, more precisely, that the 
Reynolds number of the flow is small. In the case in which the flow is two-dimensional 
or has radial symmetry, the introduction of a stream function and the Stokes assumption 
together with the no-slip boundary condition on the body and an appropriate assumption 
at infinity, yields a boundary value problem which has been solved in a number of 
instances. See Dryden, Murnaghan and Bateman [1], pp. 295-312 and the references 
of [2]. 

In the axi-symmetric case a body (or configuration of bodies) having an axis of 
symmetry is immersed in a flow which has a uniform velocity parallel to the axis of 
symmetry. It is then reasonable to suppose that the resulting flow pattern is identical 
in all planes through the axis of symmetry. If we introduce cylindrical coordinates 
x, r, 0) into the flow space, where z(— © < x < @) is measured along the axis of sym- 
metry, 7(0 < r < ©) normal to this axis, and 0(0 < @ < 27) is measured with respect 
to an arbitrary plane through the axis of symmetry, then @ plays no further role in our 
analysis because of the axial symmetry. Now let the velocity of the fluid be u(z, r) = 

u,(x, r), u,(x, r)) and introduce a stream function (2, r) by the equations 


4 
u. = fe. u, = -* =. (1.1) 


ror’ 5 


S| 


By a well-known procedure (Payne and Pell [2]; Milne-Thomson [3], pp. 521-523) the 
differential equation to be satisfied in the region of flow D is found to be 


Dy = 0, (1.2) 


££. _i2 
ox” * Or =r Or 


L., = (1.3) 
Let the trace of the boundary of the body in a meridional plane be C (Fig. 1). Then 
the condition of vanishing velocity on C can be stated in the form 


> 


y= °| (1.4) 
on C 
2 ’ 
el (1.5) 
on J 


where n is the unit normal to C exterior to the body. 
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If the uniform velocity of the flow at infinity is u = (U, 0), then y must satisfy the 
condition 
lim y = }°°U + O(p), (1.6) 
where 
p=2+r. (1.7) 


II. Representation of the solution. It is expedient to define a second stream function 
y, by the relation 


y = 3Ur’ — y, . (2.1) 
We then find that ¥, must satisfy the equation 
L,y, = 0 (2.2) 
in D, subject to the conditions 
¥, = 3Ur, (2.3) 
OY; - ’ 
oes = Ur ie? (2.4) 


on C, as well as the condition that y, give rise to a vanishing velocity at infinity. 
Following Weinstein [4], solutions of 


Ov dv , kav hai a 
L.v) = —3 + = -—=(0 (2.5) 
OX or Tr or 


k real, are known as generalized axially symmetric potentials, and denoted by y*. Payne 
has shown that in certain regions any solutions of the repeated operator equation 


[5 
(2.2) can be represented as a linear combination of any two of the functions 
. 2 2 2 »3 2.2 23 
ary, b) ary’, c) pry, 
2,1 45 
dry, e)ry. 


In a previous paper [2] the authors have used this theorem to obtain the solution of 
certain problems in Stokes flow; it will be used here to obtain still another. 

III. The flow about a spindle. We introduce bipolar coordinates (€, 7) into the 
(x, r) plane by the transformation 


z = ib cot (¢/2), (3.1) 
where z = x + ir, ¢ = + 7y, and b > Oisa constant. In terms of the (&, 7) coordinates 
—_ _bsinh 7 | = bsinté _ (3.2) 

cosh » — cosé cosh 7 — cos&’ 


where — © < 4 < » andO < é < zw. We define a spindle to be an object whose surface 
is obtained by revolving the curve § = &(0 < & < 7) about the z axis. This curve is 
the arc lying in r > 0 of the circle which passes through (+), 0) and has its center at 
(0, —b cot &). The surface of the spindle is thus § = & and the region of flow 


DisO<§E§<&,-—2© << @ (Fig. 1). It is advantageous in considering flow about 
the spindle to choose é and 7 as independent variables. 
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Fic. 1. 


Using the results of [5] we represent the solution of (2.2) in the form 
nary +(e — vy. 


It is not difficult to see [6] that any ¥*"*'(é, 7) which is even in » may be expressed in 


(3.3) 


the form 


p2net = (ec s} am ag t)°? 
y = (cosh n — cos £) 


sin [ F(a)K®” (cos £) cos anda, (3.4) 
where F(a) is an arbitrary function of a, 

(3.5) 
(known as the conal function; see Hobson [7], pp. 444-453) and K{" (A) = d"K,(A)/dX". 


Substitution of (3.4) in (3.3) permits us to write y, in the form 


K, (cos) = Pj a-1/2 (cos &) 


vy, 1Ur'(s — t)'? [ [A(a)tK S(t) + Bla)K.(t)] cos an da, (3.6) 


v0 


where 
(3.7) 


s = cosh 7 t= cosé 
and A(a), B(a) are functions to be determined in such a way that the boundary condi- 
tions (2.3-.4) are satisfied. The first of these yields 


| [A(a) toK 2 (to) + Bla)K .(to)] cos anda = (s — oo. (3.8) 
Next we observe that 0( )/dn = 0 is equivalent to 0( )/d& = 0 on = & , and make 
use of (3.8) to show that (2.4) can be written in the form 


f , 
a 0 wtide , . OK (ty) 0 ~ 

e a\to/\ a en 2 2 ¢ 
| A a) a [toK a (to)] + Bla) at, } cos an da ate (s > ae (3.9) 


The conal function K,(—?) may be defined ({7], p. 446) by 


9) 2 _ @ - 
K,(-—t) = - cosh at | [cosh u — cos ¢]™ 
T 


“ O 


fo 


cos au du. (3.10) 
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But for 0 < & < 7, (cosh 7 — cos &)~’”” satisfies the hypotheses of the Fourier integral 
theorem, and utilizing (3.10) we thus have 


—— 5 a ° _ 
(s— tt)? == | COS an if [cosh u — t]7'”* cos au iu} da 
0 


T Jo 


(3.11) 


da, 





— 91/72 [ K .(—to) cos an 


cosh am 


where t, = cos & . Thus the right hand member of (3.8) may be replaced by the integral 
(3.11). For 0 < & < m one may differentiate (3.11) with respect to ¢ , and obtain a 
representation for the right-hand member of (3.9). We are led in this way to two linear 
equations for A and B whose solution is 


1/2 


( , = -—7- ral— \K : ‘é — my ee 3 0 id . 2) 
A(a) pace ape [K.(—t) Ka (to) Ka (— ty) K a(to)], (3.12) 
, 3 ‘ »(1) r(1)/s \ > , d >(1) ‘ / ‘ 
Bla) = ——— [6K 2" (tb) Ka (—b) —Kal(—to) = (boK a (to))], (3.13) 
2 cosh ar dt, 
where 
Q = te[K2(t) — Kalle) 5 (toK 2 (to). (3.14) 


More compact expressions for A and B can be obtained, however. Suitable regrouping 
of the terms of 2 and repeated use of the fact that K, satisfies a Legendre equation 


yields 
“ (a — #)2]) = —2K,()K? (0). (3.15) 
Integration of this with respect to ¢ and division of the result by 1 — ¢ gives 
o-7*5 [ K.A)K2(1) dr, -1<& <1. (3.16) 


The bracketed quantity in A is simply the Wronskian of K,(—?t) and K,(¢t) evaluated 
at t) , which (Neumann [8], pp. 207-210) is —2 cosh ax/x(1 — #5). Thus we obtain from 
(3.12) 





1/2 1 -1 
A(a) = = ([ K (pK? (7) dr) (3.17) 


for —1 < % < 1. 
By essentially the same procedure as was used in discussing Q, it can be shown that 


; 
5 \¢ — PK (0K. (-t) — K-09 £ xeon} = —-2K2()K.(—t) (3.18) 
holds for —1 < ¢ < 1. This may be integrated with respect to ¢ from t (—1 < t < 1) 
to 1 and divided by 1 — @ in order to obtain the bracketed portion of B(a). It will be 
noted that (3.18) is not valid at = +1, since these are singular points for either K,(¢) 
or K,(—1t) and their derivatives, but the relations ([8], pp. 207 and 209) 
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lim (1 — )’K2(-2d = Pr cosh ar, 7 =1,2,-:> (3.19) 
K.(1) = 1 (3.20) 


co -adeole + OTe +O] beet} 
j=1,2,--- (8.21) 


permit the integration just mentioned to be carried out up to ¢ = 1. We find then that 


3.13) yields 





B(a) = Ata),h(a? + 1) - —= / K(—) K(x) ar} (3.22) 


4 cosh ar 


for -—1 < t% <1. 

The insertion of (3.17) and (3.22) in (3.6), and of the result in (2.1) gives the stream 
function for the Stokes flow about the spindle. 

It is easily verified that for § = 2/2(t) = 0), Eq. (3.6) with A(a) and B(a) defined 
by (3.12) and (3.13) [or by (3.17) and (3.22)] gives the correct result for the flow about 
a sphere. In this case we note that 


lim K,(t) = lim K,(—t) (3.23) 
and 

lim K"(t.) = —lim K"’(—%#). (3.24) 

te-0 te-0 


For ¢, = 0, (3.12) and (3.13) then gives 
A(a) = —2°”/cosh ar, (3.25) 


Bla) = 2'”/cosh ar. (3.26) 


From (3.6), the expression for ¥, becomes 





> rly -_ »(1) 
vs = 270r@ — 97 [BESO — ths OI aos an da. (3.27) 


im cosh ar 


The integral on the right hand side of (3.27) is easily evaluated from (3.11) and the 
expression obtained by differentiation of (3.11) with respect to & . We merely replace 
£, by + — £ in the resulting integral formulas and obtain at once the following expression 


7 


for y; : 


LUr(s — t)'7{(s +O? + ts + 7}. (3.28) 


Il 


V1 
This may be rewritten as 


= 1Ur(s — (366 + 0"? — (6 — Os + 07} 


4Ur* {3bp™* — b°p™*} (3.29) 


where b is the radius of the sphere. This is the well known result of Stokes (see [3]). 
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IV. Drag of the spindle. It was shown in [2] that the drag P of an axially symmetric 
body for which the flow region D is simply connected is given by 


—s py (4.1) 


Substituting from (1.7), (3.2), and (3.6) in (4.1) we obtain 


? fi ee poets —S 
ae = ue lim (s + #)'” | [A(a)tK?’(t) + B@)K.(t)] cos an da, 
Siu Z g 0 J0 


9! 2 a @ 
= Ub | [A(e)K (1) + Ble)K,(1)] da. 


v¥Q 


II 


From (3.20-.21) we obtain K,(1) and K;’’(1), and thus the drag of the spindle becomes 


P = 2°"eubU | [2B(a) — (a? + 4) A(a)] da (4.2) 
or 

P = SxubU [ te aa, (4.3) 
where 


F(a) = | 


With A(a) and B(a) defined by (3.25) and (3.26), Eq. (4.2) gives the well known result 
for the sphere. 

The reviewer has kindly called our attention to the fact that tables of the functions 
K,(t) are now being compiled (M. I. Zhurina and L. N. Karmazina, Tablitsky funktsit 
Lezhandra P_;/2+;,(a).) These tables should facilitate the computation of the drag coeffi- 


K.(—n Ke’ (7) ar / | K.(n)K (7) dr. (4.4) 


cient P as a function of & . 
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ON PARAMOUNT MATRICES* 


BY 
PAUL SLEPIAN 
Hughes Aircraft Co., Culver City, California 


Introduction. A paramount matrix is a symmetric matrix of real numbers such 
that any principal minor must be at least as large as the absolute value of any other 
minor of the same order built from the same rows. In this paper a new property of 
paramount matrices is cited and proved, and some remarks are made relevant to the 
problem of the synthesis of multiport resistive networks. Unfortunately, to achieve a 
concise and rigorous presentation, a rather formal mathematical notation was developed, 
and some of the results, although not particularly profound, may be obscured by the 
mathematical formalism. Hence, a long introduction is included which describes the 
results in less exact but more comprehensible terms. 

After several paragraphs introducing the notation, it is shown in Theorem 9 that 
if any elements in the main diagonal of a paramount matrix are increased and all other 
elements remain fixed, the resulting new matrix is still paramount. The remainder of 
the paper is devoted to the problem of the synthesis of multiport resistive networks. 

It is well known [1] that if A is either the open-circuit impedance matrix or short 
circuit admittance matrix of a multiport resistive network, then A is paramount. How 
ever, it is not known whether paramountcy is sufficient for the realization of such a 
network. More precisely, if A is an arbitrary** paramount matrix of order m, it is not 
known whether there exists an m-port resistive network such that A is either the open- 
circuit impedance matrix or short-circuit admittance matrix of this network. Tellegen 
and Elias [2] have shown that the preceding statement is true if m < 3, but when m is 
arbitrary, the answer to the question remains one of the leading unsolved problems. 
This paper does not attempt to answer the question for an arbitrary m, but it does 
indicate a technique suggested by the the preceding result which may ultimately help 
to solve the problem. 

Observe that an arbitrary mth order symmetric paramount matrix A has m(m + 1)/2 
independent elements. Thus, an m-port resistive network whose open-circuit impedance 
matrix or short-circuit admittance matrix is A can be expected to contain, in general, 
at least m(m + 1)/2 resistive branches. 

Note, however, that by Theorem 9, the set of all real numbers xz such that a para- 
mount matrix results when the element A,, of A is replaced by 2, is a semi-infinite closed 
interval bounded on the left. Let b be the least element of this interval and let B be the 
matrix obtained from A when 4A,, is replaced by b. Now, if there exists a resistive net- 
work whose open-circuit impedance matrix is B, it is reasonable to expect that there 
exists such a network with (m(m + 1)/2) — 1 resistive branches since the number of inde- 
pendent elements of B is (m(m + 1)/2) — 1. Furthermore, if such a network is available 
and a resistance of A,, — b ohms is introduced in series with port 1, then a resistive 


*Received December 26, 1959. 
**An arbitrary paramount matrix of order m is a paramount matrix with exactly m(m + 1)/2 


independent elements. 
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network is obtained whose open-circuit impedance matrix is A. Thus, to realize A, we 
need consider only geometrical structures realizing B, with (m(m + 1)/2) — 1 resistive 
branches. This idea can be generalized as follows. 

In this paper, in Sec. 17, a number n(A) is defined for any paramount matrix A. 
This number (A) is the maximum number such that the realization of A can be achieved 
by combining n(A) resistances in series or parallel with the ports of a resistive network of 





[order (A)][1 + order (A)] _ 44) 


2 
resistive branches that realizes a certain minimal paramount matrix obtained from A. 
This minimal paramount matrix obtained from A is minimal in the sense that it is 
irreducible as defined in Sec. 14. Thus, realization of A has been simplified to the con- 
sideration of the realization of such a minimal paramount matrix by a resistive structure 
with only 


forder (A)]| = order (A)] _ n( A) 





resistive branches. This represents a simpler geometrical realization problem, and 
because of the maximality of n(A), it is the best simplification which can be effected 
by utilizing the concept described above. 

It is easy to show that for any arbitrary third-order paramount matrix A, n(A) = 3. 
Thus as shown by Elias and Telegen [2], any arbitrary third-order paramount matrix 
A is alwuys realizable as the open-circuit impedance matrix of a resistive network of 
six branches which has three resistive branches in series or parallel with the ports and 
three more resistive branches realizing the minimal* paramount matrix obtained from 
A. The same geometrical structure will always produce a resistive three-port network 
whose open-circuit impedance matrix is equal to an arbitrary third-order paramount 
matrix A. 

Unfortunately, the situation is more complicated in the case of an arbitrary fourth- 
order paramount matrix. It is not true that for each arbitrary fourth-order paramount 
matrix A, the number n(A) is the same. It is easy to exhibit such a matrix A with n(A) =2 
and another such matrix A with n(A) = 4. Clearly the minimal paramount matrix 
resulting in the former case will have eight independent elements, while the minimal 
paramount matrix resulting in the latter case will have six independent elements. This 
leads the writer to conjecture that if paramountcy is indeed a sufficient condition for 
the realization of a fourth-order matrix as the open-circuit impedance matrix of a re- 
sistive network, the same geometrical structure will not suffice for each arbitrary fourth- 
order paramount matrix A but will depend upon the number n(A). The formal math- 


ematics follows. 


*Consider the subnetwork of Tellegen’s network which does not contain any of the three resistances 
in series or in parallel with the ports; agree to call any three-branch resistive three-port network a 
minimal resistive three-port, if and only if it has the geometrical configuration of this subnetwork. It 
was observed by Prof. R. M. Foster of the Polytechnic Institute of Brooklyn that a third-order para- 
mount matrix is realizable as the open-circuit impedance matrix of a minimal resistive three-port net- 
work, if and only if the matrix is irreducible (as defined in Sec. 14). 
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1. Definitions. 
R {x | x is a real number} 
il = {x | xis a positive integer} 
{x | for some y, (x, y) ef} 
iv) rng f = {x | for some y, (y, x) ef} 
f{A=f\ {(z, y)| ze A} 
If x e w, then 
(x) = wl \fyly <2} 
pA the cardinal number of A. 


2. Definitions. 
For each x ew 
S* = {u | uw is an increasing sequence with 1 < pu and rng u C (2)}, 


For each x e w and 7 e (x) 


Si = SC) fultcgmg u}. 


3. Definitions. 
VW = {A|A isa function, dmn A = ((m) X (m)) for some mew andrng A C R}. 


If 4 e M and (7,7) edmn A, 


i 


1 
Ll 


iii) If 4 eM, order (A) = (p dmn A)'” 


iv! K = M(\ {A]|A,, = A,; for (¢,7) e dmn A and 2 < order (A)}. 


For m €, 


M”" = M(\ {A|m = order (A)}. 


If A e M, det (A) is the customary determinant of the square matrix A. 
vii) If A e M and det (A) ¥ 0, then A“ is the customary inverse of A. 


4. Definitions. 
fu * x] 
iat increasing sequence of positive integers such that 


dmn[u « x] = (x — pu) 


inf ((x) — rng u) 


lu + 2) 
[fu « x]; = inf ((x) — rng u — rng([u * x]|{i — 1))) 
forte (x — pu) and i #1. 
5. Definitions. 
For ieM 
G(A) = (S7e"™ x Sortet“*”) CV {(u, v) pu = po}. 
For A e M and 7 e (order (A)) 
GA) = G(A) 1 (Su K Serer), 
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6. Definitions. 


For A e M and (u, v) e G(A) 


A“ is that element of Jfo™der‘“’~™ 


such that 


A*),; = A({u * order (.4) ], , [v * order (A)],) 
for 


i, j) e (order (A) — pu) X (order (A) — pu)). 


7. Definition. 


A is paramount if and only if A e A and 


det (A%) — [det (A‘)} > 0 
whenever (wu, v) e G(A 
8. Definition. 


If Ae M,2eRand7e (order (A 


is that element of M/erder'*’ such that 


if (j, k) edmn A and (j, k) 


9. Theorem. 


Suppose that A is paramount and that 7 e (order (A)). Let 


B = {al 


x|[A, 7, x] is paramount}. 
Then B is a closed, semi-infinite interval, bounded on the left. 
Proof. 


Let 0 < te R. It suffices to show that [A, 7, A;; + Jt] is paramount. Pick any 
(u, v) ¢ G,([A, 7, Ags + 7). 
It suffices to show that 
det ([A,7i, A,; + t]3 > [det ({A, 7, Ais + tf]. 


Suppose first that 


pi = p= order (A) — |. 
In this case 


det ([A, 27, Ai + th) = Au + 8. 
Also, if 7 g rng v, then 
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det ({A, t, A;; + t]:) - Aj; + t; 
while if 7 e rng v, then 
det ({[A, 7, Ag + i) = Ai; 


for some j e (order (A)) with j ¥ 7, and paramountcy of A implies the desired inequality. 
Thus, assume now that 


pu = po < order (A) — 1. 
Since we Sot" 7 e rng [uv * order (A)]. Let 7 e dmn [wu « order (A)]. such that 
fu « order (A)]z = 7. 


Suppose now that 7 g rng v. 
Let 7 e dmn [v « order (A)] such that [v « order (A)]; = ¢. Then 


det ((A, 7, Ay; + ) = —I*itdet ({A, i, Ass + O)\3}) + det (4%). 
Clearly there exists (a, 3) e G(A) such that 
({A, 4, Ass + &Dia’7; = Ai. 


Thus, 


det ([A, i, Ai; + tJ) = —1°*7tdet (A®) + det (A‘). 


Note that the above equation is valid when u = v andi = j. Hence, the fact that A 


is paramount implies 
det ({[A, 7, Ai; + td) — [det ([A, 7, Ais + é]D 
> tdet (A®) + det (A*) — (¢ |det (AD| + det (A%)|) 
= (det (A) — /det (A®)}) + (det (A%) — |det (A%)]) > 0. 


Finally, if 7 e rng v, then 
det ({A, i, Ai: + tJ) = det (A¥), 
and 
det ({A, i, Ay; + ti) — |det ({[A, 7, Ai; + 8! 
= tdet (A%) + (det (A%) — det (A%))) 


IV 
o 


The proof is complete. 
10. Definition. 
B is a reduction of A if and only if A is paramount, and 
B = [A, i, inf {t |[A, 7, ¢] is paramount} } 
for some 7 ¢ (order (A)) such that 


A;,; > inf {t|[A, 7, ¢] is paramount}. 
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11. Definition. 
P is a reduction sequence of A if and only if A is paramount, and P is such a sequence 


that 
and 
P, is a reduction of P;_, for 1 < 7edmnP. 
12. Remark. 
It is clear from Definitions 10 and 11 that if P is a reduction sequence of A, then 
pP < order (A). 
13. Remark. 


Suppose that 7’ is such a sequence that 
Foriedmn TJ T(z) is a sequence to M. 


(i) 
(ii) 7(1) is a reduction sequence of A. 
(iii) (7'(z)),rq) is non singular for each 7 e dmn 7 
iv) T(¢ + 1) is a reduction sequence of [(7'(7z)),7:;,)) ° for each 7 e dmn JT such 
| 4 AC) ) prea 
7+ 2)edmn 7, 


iv) 
that (¢ + 1) edmn T. 
Then, clearly, by Jacobi’s theorem, whenever 7 e dmn T and (7 + 2) 

pT(t + 2) < pT(i) < order (A). 
Thus a sequence 7 satisfying the above conditions must be finite. 


14. Definition. 
A is irreducible if and only if A is paramount and there exists no B such that B isa 
reduction of A. 
15. Definition. 
T reduces A completely if and only if A is paramount and T is such a finite sequence 
that 
(i) ForzedmnT T(z) is a sequence to M. 
(ii) 7(1) is a reduction sequence of A. 
(iii) If 1 < pT then (7(2)),7;,) is non singular for 7 e (pT — 1) 
(iv) If 1 < pT then T(i + 1) is a reduction sequence of [(7(z for 
ie (pT — 1). 
(v) (T(pT),r~@r) is irreducible and singular, or, 1 < p7', and (T(pT)),7,.,r, and 
(T(pT — 1)),rq@r-1y are both irreducible. 
16. Remark. 
If T reduces A completely, then 7 is maximal in the following precise sense. 
Theorem. 
Let T reduce A completely and let U reduce A completely such that 7 C UU. Then 
T = U. 


17. Definitions. 
If T reduces A completely, 


(i) , 
nT, A)= > pT). 
iedinnT 


to 
o 
© 
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(ii) If A is paramount and A is not irreducible, 
n(A) = sup {n(T, A) | T reduces A completely}. 


If A is paramount and irreducible n(A) = 0. 


18. Remark. 
It is not known whether the following statement is true: 
if T reduces A completely and U reduces A completely, then 


n(T, A) = n(U, A). 


19. Remark. 

If A is paramount, the number n(A) defined in 17 (77) may be relevant to the problem 
of the realization of A as either the open-circuit impedance matrix or the short-circuit 
admittance matrix of a resistive network. This topic is discussed in the Introduction. 
Clearly, if 18 is true, the computation of n(A) will be greatly simplified, and thus de- 
termination of the validity of 18 may be of help in the ultimate solution of the realization 


problem. 
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BOOK REVIEWS 


(Continued from p. 244) 


of a competitive economy, in particular, the convergence over time of the price vector to an equilibrium. 
For the latter problem a key condition is gross substitutability of commodities, which in mathematical 
terms says that a certain matrix is of Metzler type (i.e., all nondiagonal elements are nonnegative. ) 

Volume II is devoted to continuous games, in particular polynomial games, games with convex 
or bell-shaped kernels, games of timing (duels) and poker models. The main emphasis is on finding the 
specific structure of the optimal strategies for each player. The art of doing this is somewhat like that 
of solving differential equations. In each case one must judiciously choose a method which exploits 
the special features of the kernel; and quite often a good initial guess about the form of the solution 
is also needed. The methods involve, for example, the geometry of moment spaces, integral equations 
with positive kernels, and the Neyman-Pearson lemma. 

These books are intended partly as textbooks for graduate students and mature undergraduates. 
They ought to have an important influence on the coming generation of economists and operations 
researchers. Numerous good problems have been included. In addition the reader is sometimes expected 
to fill in not entirely obvious details as he goes along, e.g., (p. 123, Vol. I) the fact that a linear function 
which is bounded above on a closed, unbounded, convex polyhedron P attains its maximum in P. There 
are Appendixes (worth reading for their own sakes) dealing with vector spaces, matrices, convex sets, 
and convex functions, to fill possible gaps in the reader’s background in these subjects. A knowledge of 
the rudiments of general topology and real function theory will prove helpful, especially in Volume II. 
The two volumes are independent of each other; in fact, both have the same first chapter (about matrix 


games) and the same Appendixes. 
WENDELL H. FLEMING 


Hydrodynamics. By D. H. Wilson. St. Martin’s Press, New York, 1959. viii + 149 pp. 
$5.50. 


The purpose of this book is best described by quoting a paragraph from its preface. 

‘This book is intended primarily for students specializing in mathematics or theoretical physics 
but will also be found useful by General Degree and Engineering students. Its object is to provide a 
brief and concise introduction to ‘classical’ hydrodynamics while at the same time providing the basis 
and acting as a starting point for the study of the many modern developments in fluid mechanics.’ 
The author suceeds admirably in reaching his objectives. He uses vector calculus to derive, in a 
1 manner, the equations of fluid mechanics in Chapter I. His representation of the basic princi- 


very lucid 
be recommended for perusal to students who value a fortunate combination of physics and 


ples can 
rigorous mathematics. Chapter II proceeds, in conventional fashion, with a study of two dimensional 
motion of the inviscid incompressible fluid. Chapter III consists of a careful study of vortex motion. 
Chapters IV and V are devoted to the conformal mapping technique for two dimensional potential 
flows and the stream function technique for axisymmetrical potential flows respectively. The book 
closes with a study of viscous motion in Chapter VI. It contains 147 pages. 

The author is to be commended for this clear and concise presentation. The book is a worthwile 


addition to the library of any student of fluid mechanics. 
P. F. MAEDER 


to 
“J 
— 


THE UNIQUENESS QUESTION FOR WAVES 
AGAINST AN OVERHANGING CLIFF* 


BY 
R. SHERMAN LEHMAN 


1. Introduction. For the problem of three-dimensional waves over sloping beaches 
Roseau [11] and Peters [10] have found explicit solutions for all slope angles. These 
solutions are in the form of integrals in the complex plane resembling those obtained in 
the inversion of the Laplace transform. It is not, however, clear that all solutions satis- 
fying appropriate conditions can be represented in this way, and indeed in general the 
uniqueness question is still open. The purpose of the present paper is to study the unique- 
ness question in a particular case. 

The problem of surface waves over sloping beaches when treated by the linearized 
theory leads to the question of determining a velocity potential @(x, y, z, ¢) which is a 
solution of Laplace’s equation in the space variables z, y, and z and satisfies two different 
boundary conditions on different parts of the boundary. Suppose the z-axis is taken 
along the shore, the y-axis is directed vertically upward with the free surface at y = 0, 
and the x-axis is directed outward from the shore. Then on the free surface y = 0, x > 0 

a 4. Ob _ 0 

at Ty 
where g is the acceleration due to gravity. On the bottom the normal derivative 
d@/dn = 0. For a uniformly sloping beach with angle 7a between the surface and bottom, 
the bottom is given by the equation @ = — ma, when z = rcos #,y = rsin@,r > 0. 
If a = } the “bottom” becomes a vertical cliff, while for } < a < 1 the “bottom” becomes 
an overhanging cliff. Our considerations here will be primarily concerned with the case 
a 3/4. 

Restricting consideration to motion which is simple harmonic in the time ¢ and the 
space variable z, one seeks solutions of the form 


&(x, y,z, 2) = ¢(x, y) cos kz cos (ot + 8). 


Progressive waves having crests which at infinity make an arbitrary angle with the 
shore line can be constructed by taking linear combinations of standing waves of this 
type. 
With an appropriate choice of units of length and time so that o°/g = 1 one obtains 
the following mixed boundary value problem: Find the solutions of the equation 

Ag(x, y) — k’g(x, y) = 0 (1.1) 


in the sector — za < 6 < 0, satisfying 


dg ‘ 
_ «2? for 6=0, (1.2) 


*Received November 12, 1959. The author wishes to acknowledge support by the Office of Naval 
tesearch under Contract Nonr-222(62). 
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se = 0 for 06 = —7a. (1.3) 

Roseau [11] and Peters [10] have found solutions for this problem for k < 1 and all 
angles ta < 7. Roseau [12] has also studied solutions for angles ra < 7/2 when k > 1. 
Some of these solutions, which include the Stokes edge waves, have been considered also 
by Ursell [15]. The question of determining all solutions of the boundary value problem 
satisfying appropriate conditions at infinity and at the origin has been studied by Stoker 
[13] and Weinstein [17] for0 < k < 1 in the case a = 1/2. (See also [14, p. 84].) A proof 
of uniqueness for a = 1/4 can be constructed by using Roseau’s special investigation 
of that case [11, pp. 27-30] together with Stoker’s and Weinstein’s methods. In this 
paper we shall use similar methods to handle the case a = 3/4. 

For i; = 0 the waves are two-dimensional. In this case solutions for all angles were 
found earlier by Isaacson [4]. Stoker [13] has treated the uniqueness question for a = 1/2n. 
3rillouét [1] has studied the uniqueness question for a = p/2q following the line of Lewy’s 
original investigation [6]. A uniqueness proof for all a < 1 for k = 0 can be constructed 
using the methods and results employed by Lewy in [8]. 

The present investigation was begun with the aim of finding out the true situation 
with regard to uniqueness for k > 0. There is no possibility of answering the question 
for all angles by the method used here, and even the case of all angles rp/2q, which in 
principle probably could be handled by the methods employed here, seems uninvitingly 
complex. One can hope, however, that knowledge of a few special cases may help in 
guiding a general consideration using methods similar to those employed by Lewy 
[7, 8, 9] and the author [5]. 

2. Reduction to a simpler boundary value problem. We wish to determine for 
a = 3,4 all solutions of the boundary value problem (1.1)—(1.3) which have continuous 
second derivatives in the sector — 37/4 < @ < 0,r > 0 and satisfy two additional 
boundedness conditions. At infinity we shall require that 


p+] 22) +] 22] <a fo r>R,, (2.1) 
| | 


where ./ and R, are some positive constants. Near the origin we shall assume that 


| 


dg | | ag | 1 —7/3+6 
nite) <& (2.2) 


2 


le| + 





< 


where C and ¢ are positive constants. Observe that (2.2) allows ¢ to have a logarithmic 


singularity at the origin but excludes solutions which behave like r-*”” at the origin. 
Because the normal derivative of g vanishes on 0 = — 32/4, the function ¢ can be 
extended by reflection to the sector — 37/2 < @ < O. In the larger sector ¢ will again 
satisfy the conditions (2.1) and (2.2). Similarly, since dg/dy — gy vanishes on 6 = 0, 
it can be extended by reflection to the sector — 37/4 < 6 < 32/2. The resulting func- 
tion d¢/dy — ¢g is multiple-valued when considered as a function of x and y, but it can 


be considered as a single-valued function of r and @. It satisfies the following two relations 


uniformly in @: 
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‘fe ' : 
?— y= Ol) for r— @, —3r 2 0< 34 
OY (2.3) 
uo 7. 
= — y= O(r-"*4 for r—-0, <0 < dr. 
Cc 
Since the function ¢ is symmetric with respect to the line 6 = — 32/4 and dg/dy — 
¢ = Ofor 6 = 0, we conclude that d¢/dx — ¢ = 0 for 6 = — 37/2. Define 
soon hi le ten - 


The function y is then a solution of the differential equation Ay — k*°y = 0 for —3x/2 < 
@ < 3x/2 because the differential operators involved have constant coefficients. Also 
Y = Ofor 6 = Oand for 6 = — 3x/2. Thus y isa solution of a much simpler boundary 
value problem. 

The method we employ is first to solve this boundary value problem for y and then 
use (2.4) to determine ¢ by the solution of ordinary differential equations. Finally we 
must verify that the ¢ obtained is a solution of (1.1)—(1.3) and satisfies the boundedness 
conditions (2.1) and (2.2). 

It is possible for all a = p/2q with p odd to perform a similar reduction to a problem 
with vanishing boundary values for a sector with angle 2ra. This simpler boundary 
value problem can be solved, (see Brillouét [1]) but since for larger p and gq the differential 
operator used in the reduction is more complicated, the determination of ¢ from the 
solution of the simpler boundary value problem presents difficulties. 

We first use (2.1) and (2.2) together with the following lemma to derive bounds 
for y at the origin and at infinity. 

Lemma 1. Let R, be an arbitrary positive number. Let u be any solution of Au — k’u = 0 
existing in a circle of radius d < R, about (x , yo) and let M, be the supremum of u inside 
this circle. Then 


M 
Ly) 








u Ou 
ae (20 5 Yo) | + | By (% » Yo) 


where L is a constant dependent only upon the differential equation and the number R, . 

This lemma is a special case of known theorems. In particular it is a special case of 
results of Gevrey [3, p. 148] and it also follows from the Schauder interior estimates 
as stated by Douglis and Nirenberg [2]. 

We now use the lemma together with (2.3) to obtain estimates for y at infinity and 
at the origin in the sector — 37/4 < @ < 32/4. For the estimate at infinity we use a 
circle of radius 1 about the point (z, y). For the estimate at the origin we use a circle of 
radius r/2 about the point (xz, y). We find that the following estimates hold uniformly 
In @: 


V(x, y) O(1) for r—> &, = 0 < 3x, 
W(x, y) = Or "”**) for r—0, —ir <0 


In the following discussion it will sometimes be convenient to consider y as a function 
of r and 6. To avoid possible confusion we set 


Vir, 0) =Y(a, y). 


rs 
m 


II 


(2.5) 


lA 


3y 
40 
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Since ¢ is symmetric with respect to the line 6 = — 37/4, the function y in view of 
(2.4) also must be symmetric with respect to @ = — 32/4. In other words 
V(r, 6) = U(r, — $n — 8). (2.6) 


By reflection across the line 6 = 0 we obtain 
Wr, d=) = —WV(r, — 8). (2.7) 


Combined with (2.6) this yields the information that WV is a periodic function of @ with 
period 37. For any fixed value of r the function Y can therefore be expanded in a con- 
vergent Fourier series which because of (2.6) and (2.7) must have the form 


V(r, 0) = >> c,(r) sin 3(2n + 18, (2.8) 


where the coefficients c,,(r) are given by 
8 0 
c,(r) = =| U(r, 6) sin 3(2n + 1)6 dé. (2.9) 
3m —34r/4 


The function V(r, @) satisfies the partial differential equation 


Differentiating (2.9) and using this differential equation, we obtain 
Ds 2, J2 4 
ty if ) a 4 - (- \ 
Cn’(r) + = eal) [ + {2 (2n + } k. r) 
i ete ‘ | wr 
ae | | 2 + (§(2n + 1))°W | sin 3(2n + 16 dé. 


The right side of this equation is equal to zero, as can be seen by integrating the first 


term of the integrand by parts twice and using the fact that ¥ = 0 for 6 = 0 
and dV/a6 = Ofor @ = — 3/4. It follows that the function c,(r) is a linear combination 
of modified Bessel functions of order 2(2n + 1)/3 that is 

Cal) = Aglocensiysa(kr) + ByKocone1/a(Fr), (2.10) 


where A, and B, are real constants. 
Next we use the estimates (2.5) to prove that all of these constants except B, and B, 


are equal to zero. For p > © over positive values the following asymptotic relations hold 
(see [16, p. 202]): 
I,(p) ~ (2mp)-"’e’, 
(2.11) 


\'2 
Kp) ~ (=) e. 
Inserting the first estimate of (2.5) into (2.9) and letting r — ©, we conclude that 
A, = Oforn = 0, 1, 2, --- . Also for p — 0 over positive values 
K,(p) = O(p"’), (2.12) 


when » > 0 [16, p.77-78]. Hence the second estimate of (2.5) when combined with 
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(2.9) allows us to conclude that B, = 0 for n = 2, 3, 4, --- . Thus we find 
V(x, y) = B.K2,(kr) sin 36 + B,K,.(kr) sin 20, (2.13) 


where B, and B, are real constants. 

3. Solution of ordinary differential equations. In Sec. 2 we proved that if g(x, y) 
is a solution of the boundary value problem (1.1)-(1.3) satisfying the boundedness 
conditions (2.1) and (2.2), then the corresponding function y (x, y) determined by 
(2.4) must have the form (2.13). In this section we shall prove that corresponding to 
each function y of the form (2.13) there is at most one function ¢ satisfying (1.1)-(1.3), 
(2.1) and (2.2), thus establishing that for each k > 0 there is at most a two-parameter 
family of solutions of our problem. The method employed in this section is essentially 
integrating ordinary differential equations. 

We shall need the following lemma. 


Lemma 2. Suppose u(x, y) is a solution of the differential equation Au — k* u = 0 
in the stripa <2 <b, — © <y < o&. Suppose also that there is a constant C, such that 


Ou 


du | 
2 ay 


| r | on? | + 





<C; (3.1) 








fory > yo,a< 2 <b. Then the function 


y 


U(z, y) = e | e~ ‘u(x, t) dt 


+@ 
is a solution of the same differential equation in the same strip. 


Proof. By differentiating we obtain 


v 2 \ 

AU — kU = e | ere] Pula, 0 + (1 — F)u(z, 0 dt + u(x, y) + yt y), 

as is easily justified in view of (3.1). Transforming the integrand by using the differ- 
ential equation and then integrating by parts, we obtain 


. a -— 3° a 
AU —kU =e’ [ e (. — =) dt + u(x, y) + By Ue y) 


J @ 


v “ -:(Ou du am 
, [ie (% at) at = 0. 


Il 


The inequality (3.1) insures that the boundary terms at the limit + © vanish. 

Now let us make a cut along the positive y-axis of the z, y plane and restrict @ to the 
range — 37/2 < @ < 2/2. In this way for given values of B, and B, in (2.13) we obtain 
a single-valued function y(z, y) defined in the entire cut plane. Let 


v(x, y) = (2 - 1 ole, y) (3.2) 


so that v is a solution of the equation 


(2 — 1)(e, ) = vor, 0. (3.3) 
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Integrating we find 


y 
v(x, y) = e’g9(x) + e” [ e ‘Y(a, t) dt = e’g(x) + v*(a, y), (3.4) 
where g is a function of x alone. By Lemma 2 v* (x, y) is a solution of the differential 
equation (1.1) for-— © <2<0,— » <y< = andalsoitisasolution for0 < z < ~, 
— o <y< o, As we shall see, »*(z, y) is discontinuous on the entire y-axis. Since the 
function v is also a solution of (1.1), e’g(x) must be a solution, 1.e. 


g(x) + (1 — k’)g(x) = 0 (3.5) 
for the intervals — ~ < x < Oand0 < x < o. We shall denote by g,(x) and g.(x) 


the analytic functions of x coinciding with g for x < 0 and x > 0, respectively. 
By (2.11) ¥(z, #) is exponentially small for «* + ¢° large. Consequently for zx > + 


with y fixed 

v(z, y) = e’g(xz) + of 
Also by (3.2) and (2.1) we know v(z, y) = O(1) forz > + ~, y < 0. When k < 1, this 
does not impose any additional niles on g2(x); but if k > 1 among the solutions of 


(3.5), a, e’* + a, e*** with k, = (k* — 1)”, we must reject all those for which a, ~ 0. 
Siues we obtain the additional condition 

93(0) = —ky,g(0), k, = (kh? — 1)”, (3.6) 
which must be satisfied in order to obtain a ¢ satisfying (2.1). Similarly for k = 1 we 
must have g3(0) = 0. Thus (3.6) must be imposed for k > 1. 


Next, solving (3.2) for g(x, y) in the lower half plane, we obtain 
g(x, y) = &h(y) +e [ e- ‘v(t, y) dt (3.7) 


By an argument using Lemma 2 with x and y interchanged one can show that the second 
term on the right hand side is a solution of the differential equation (1.1) for y < 0. 
It follows that for y < 0 

h’(y) + - k*)h(y) = (), 


By (1.2) we must have dg/dy — ¢ = O for y = 0, x > O. For all y < 0 we have 


FY : az P 
s —g=e(h'(y) — hy) +é& | e “(2 _ 1 oct y) dt 
) 7 (3.8) 


Il 


E(') — Wy) +e | e'W(t, ») dt 


Keeping z fixed > 0 and taking the limit for y — 0°, we conclude that we must have 
h’(0) — h(O) = 

O(1) fort — + ©, we obtain from (3.7) 

g(x, y) = eh(y) + O(1) 

forz— + o with y fixed. To have ¢ bounded as required by (2.1) we must have h’(0) = 

h(0) = 0, and hence by the differential equation h(y) = 0 for all y. 


Since, as we have seen, v(t, y) = 
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The boundary condition (1.3) also must be satisfied. For y = x, « < 0 we have, in 
view of (3.2) and (3.8), 


9299 _ 99 _ Ip 


on Ox 0 y 


= [¢(x, y) + o(z, y)] — E [ e ‘W(t, y) dt + o(z, »)| 
sa (3.9) 


= (2, y) — é | e ‘yy, t)dt (since y(t, y) = Wy, 0) 


= g(zje", by (8.4). 


Consequently we must have g,(x) = 0 for all z. 
Since v and év/dx must be continuous across the negative y-axis, we must have for 
0 


q(0) = lim | e'y(x, t) dt — lim [ e ‘V(x, t) dt, 


sit: sitesi (3.10) 
9(0) = Lim [ “ma~ lim [ on OY ay. 
z-0- « Ox z-0* t-) Ox 


The function g(x) is uniquely determined by g,(0) and g3(0) since it is a solution of the 
differential equation (3.5). In the next section we shall establish that these limits exist 
for every function y of the form (2.13) and are independent of y. Without further work, 
however, we can conclude that if there is a solution corresponding to y, it is given by 


the formula 


g(r, y) = al | e ‘g(t) dt + &*" [ / e "Wr, t) dt dr, (3.11) 


< 


where 


g(x) if «>O0. 


Observe that by (3.10) and (3.11) ¢ is uniquely determined by y and also that the 
dependence is linear. From (2.13) we see that for each k > O there are at most two 
linearly independent solutions of the boundary value problem satisfying the conditions 
(2.1) and (2.2). In fact, for k < 1 Roseau [11] and Peters [10] have found two linearly 
independent solutions, one of which is finite at the origin, and the other of which has a 
logarithmic singularity there. If one used their work, to complete the proof that for 
each y the formula (3.11) represents a solution of our problem, one would only need to 
verify that their solutions satisfy the conditions (2.1) and (2.2). However, we shall not 
do this, but instead we shall give a proof which is independent of their work. At the 
same time we shall complete the study for k > 1, where it turns out that (because of 
(3.6)) there is only one linearly independent solution of the problem and it has a log- 
arithmic singularity at the origin. 

To summarize, we have proved that for each y of the form (2.13), if the limits (3.10) 
exist, then the function g(x, y) given by (3.11) is a solution of the boundary value 
problem (1.1)-(1.3). In the next section we show that the limits (3.10) exist and at the 


ih oe E if «<0 
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same time evaluate them explicitly by using Bessel function formulas. Then in the last 
section we study the behavior at infinity and at the origin of y(z, y) as given by (3.11) 
and in doing so verify that it satisfies the conditions (2.1) and (2.2). 

4, Evaluation of some limits. We shall use the formula (see [16 p. 388]) 


7 sinhv@ 
sin yr sinh B 





[ et" *K,() dt = 


#0 


valid for 0 < »v < 1, Re cosh 8 > —1, —7/2 < Imcosh8 < 72/2. It follows that 





T sinh v8 (4.1) 


| e K,(kt) dt = i sin yr sinh B 


for0 <<» <1,0<k < ,cosh8 = 1/k withO < Im cosh 8 < 2/2. More precisely, 
if k = 1, then 8 = 0; and the ratio sinh »@/sinh 8 is defined to have the value v. For 
k > 1, 6 is complex and therefore sometimes it may be more convenient to use the 
alternative expression 
vy sin vy 
k sin vr sin y 
with cosy = 1/k,0 < y < 2/2. 
To abbreviate, we introduce a limit operator 6 defined by the equation 
of = lim f(z, y) — lim f(z, y) 
z—0*T z~0- 
for any function f for which the limit exists. Then (3.10) becomes 


aw 


6 | e ‘y(zx, t) dt, 


g2(0) 
(4.2) 

7 ee | 

g0) = 6] < yx, t) dt 


“vy 


provided the limits exist. We shall prove that these limits do exist for every y < 0 and 
are independent of y. Throughout this section we shall assume that y < 0. 

We investigate separately the contributions due to each term of y as given by (2.13). 
In each case in the following discussion where we obtain a limit of an integral by passing 
to the limit under the integral sign, this can be justified by the Lebesgue dominated 
convergence theorem. First we have 


6 fe 'Krp(kl? + 2°)') sin 90 dt 
vs 


oe sinh (26/3) 
- [ e~'Kaa(kt) [sin (/3) — sin (—x)] dt = FEL, 


Using the formula 


is) 


7 k 


K,(kr) = 2 K.(kr) sin 26 


Ox OY 


and integrating by parts, we find 
L, =6 | e'K,(k{? + x*]'”) sin 20 dt 


vy 
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a ae 2 271/2 
oe s eee d 
sf a at ax MotAlt wl) a 


— “2 =-§ + ff 42 2451/2 
7 3 | we KSkle + 2°)") dt. 
Also (see [16, pp. 79, 80]) 
0 k x 
55 Kolker) = —— K,(kr) = —3 + 2p0), 


where there is a constant C such that for0 << r< © 


| pr) | < C(L + | logr )). 
Consequently 


L=-6f Z e's siaatol Zz ae"'mp(( + 2°]")dt. 


The second term is 0 since the limit of the integrand is 0 for z — 0. The first term can 
be evaluated by introducing r = ¢/z as a new variable of integration. Since y < 0, we 
obtain 





“ n x 4 is res 
4 = — —- : a eee | d =—_ — TS 
/ é [ k° ? - + t . (° bare w/z 1 + r dr 
4r 
= Re 





From the formula (see [16, p. 79]) 
zK;y(z) + vK,(@) = —2K,-,(2) 
one can derive the following identity: 


Ps ‘ 0 . 
= (K,(kr) sin v6) = ay (K,(kr) cos v6) — kK,_,(kr) sin @ — Dé. 


Using this identity and the fact that K_,(z) = K,(z), we obtain 


ine ‘ 
5 | "= {Ka(k[ + x°]'”) sin36} dt 


= 8 f ot 2 (Kip(kit + 21) costo} at 
+3 / ke'Ky (k(t + 2°]'”) sino dt = L, + Ly. 


The limit Z, can be evaluated by first integrating by parts and then passing to the 
limit under the integral sign to obtain by (4.1) 








L, = [ e'K2,3(kt)[cos (x/3) — cos (—7)] dt = 
Also 


= kf &'Kialk VE + 2)lsin (4/6) — sin (—7/2)] dt = 3°74 SE). 
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We have 
dig ‘ re] » 9 2741/2 . 
L, = 8 | ee an (K.(k[ + 2°]'’*) sin 26) dt 
ai 


vy 


= 6 - 2e7" ae. es .f42 nee 11/2 
= 0 | k? at ax (K(k t ~- 4 ] )) dt 


Ces eee an 
= § | 2 e- ay? K,(k[t' + 2” ]'”*) dt. 
Since K,(kr) is a solution of AK, — k’K, = 0, we obtain 
a 27 3” 


~~ Fe P42 ,241/2 
EF pp Kolklt? + 2°]") dt 


Ly= 6 [ 2K (kl? + 2°]) dt — 8 


= § [ (2 — 2 eK (kl + 2”]'”’) dt = 0, 


since the limit of the integrand is the same for x — 0° as it isforrt-— 0°. 
Thus 


= -B.= 
sinh 8 ‘ke sini 
(4.3) 


-B, ‘s * sinh (28/8) , 91/2_ sinh 90a) 


7 sinh (28/3) 4m 
“t 


g2(0) = B 








’(0) : : ‘ = 
92M, k sinh B = sinh 8 


where cosh 8 = 1/k, Re 6 > 0,0 < Im8 < 7/2. It is easily verified that the coefficients 
of B, in these expressions never vanish. 

When combined with (3.6) these formulas establish that for k > 1, By is determined 
uniquely if B, is given. For example for k = 1 we must have B, = 0. Thus for k > 1 
there is at most one linearly independent solution satisfying the conditions (2.1) and 
(2.2). 

5. Behavior at infinity and at the origin. We now investigate the behavior at infinity 
and at the origin of g(x, y) as given by (3.11), (2.13), (4.3), and (3.5). We begin by 
estimating 


v*(z,¥) = vi e ‘yY(zx, t) dt (5.1) 
(see (3.4)). If y is any one of the functions given by (2.13) and p is an arbitrary positive 
number, then there is by (2.11) a corresponding number c, such that 
| Wiz, y) | <ece for r= (22+ yy)” > p. (5.2) 
(It is necessary to omit a neighborhood of the origin because of the singularity of K, 
there). For y > 0, 2° + y’ > p’ we obtain 





ly*(x, y) | <e i a ass ray tae (5.3) 
Fory <0,2° + yy’ > p’, x ¥ 0 we write 
v*(z, y) = e” [ e ‘¥(zx, t) dt + e” | e‘Wr,)dt=T7,4+T.. (5.4) 
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By the results of the last section we know the integral in T, approaches a limit as z — 0° 
and also approaches a limit as z + 0°. Furthermore by (5.2) the integral in 7’, approaches 
Oas|a|— o. Hence there is a number c, such that | 7, | < ce”. 

lor the second term 7’, we have by (5.2) 


v 
| T. | < ce" [ € ‘dt|, 
J-p | 





which is less than or equal to 


c,e” 


[k-1| 
for k * 1 and is less than or equal to 
ce(| y| + p) 
for k = 1. In either case, if we define 


(Ke for k<1 


or + —— 


w= )jl-— 7 for k=1 
, for k> 1 


with » an arbitrarily small positive number, then there is a constant c; independent 
of x and y such that | T, | < c,e””. Combining this estimate with that for | 7, | and (5.3) 
we find 

lv*(z,y)| <ece"'' for +¥0, 2+y'>p’. (5.5) 


Since for y < 0, 


g(z, y) = ad | e‘g(t)dt+é [ e ‘v*(t, y) dt, 


“ 


we obtain immediately 





.? | 
g(x,y) — id e'g(t) dt | < ce”. 


| 
For k < 1 with k, = (1 — k’)'” we have 
g(t) = g2(0) cos kat + kz'g2(0) sin ket, 


and hence for x > 0 


é a e‘g(t) dt = 7 


” 
J. @— Pk) {—(g2(0) + g2(0)) cos kx (5.7) 


+ (kig2(0) — gi(0))kz' sin ker}, 


while for x < 0 since g(t) = 0 fort < 0, 
az 0 = 
| cgi dt=e / e~'g,(t) dt = eo FF (g2(0) + gi(0)). (5.8) 

If k > 1, assuming (3.6) is satisfied, we have 


kit 
’ 


g(t) = g.(Oe7 











282 R. SHERMAN LEHMAN [Vol. XVIII, No. 3 


where k, = (k’ — 1)‘. Hence for x > 0 


az —- _ — g.(O0)e~*** ‘ 
€ I e g(t) dt = sk (5.9) 
while for x < 0 
z ; i _ \ 2? 
e [ e g(t di = —92(0e . (5.10) 
Jo ] —_— ky 


Thus we see that for y < 0, ¢(z, y) is bounded outside a neighborhood of the origin 


for every function y given by (2.13) in case k < 1; and the same is true for k > 1 provided 


(3.6) holds. 
To study the behavior on the free surface y = 0 it is desirable to use a better estimate 


for v* there. By (5.2) for z > p we have | y(z, f)| < c,e**. It follows that | v*(z,0)| < 
c,e ** and hence for y = 0 the right side of (5.6) can be replaced by c,e-“*/(k + 1). 
Equation (5.7) then shows that for k < 1, g(x, 0) approaches a sinusoidal wave 
aszx—-+ o. Similarly fork = 1, ¢(z, 0) approaches a constant as x — + , the con- 


stant being — g,(0). Fork > 1, g(x, 0) = O(e"**) ast 9 + @. 
To complete the verification that condition (2.1) is satisfied we must 
| dg/dx | and | d¢/dy | also are bounded. By (3.2), (3.4), and (3.8) we have 


verily that 


00 , — ‘ , ae 
a = 9(2, y) + o(z, y) = o(z, y) + e’g(z) + v*(z, y), 


& 


Ox ‘“ 
§.11) 
00 , : Se aga” F Aap 
> = 9%, y) te | e ¥(t, y) dt = g(x, y) + v*(y, 2) 
oY J@ 
since ¥(y, ¢) = W(t, y). The boundedness then follows in view of the estimate (5.5) for v*. 
Next we study the behavior of ¢ near the origin. In view of (2.13), it is convenient 


to consider separately two different functions y, 


Yor, ¥) = Kz,3(kr) sin (26/3), 


is) 


be 


, Z ; a — 
Vi(z, y) = = K.(kr) sin 26 = a K,(kr), 
y 


Ox 


bo| 


and obtain estimates for the corresponding functions v% and v%. 
We have 


v4(z, y) = é | e'K2,,(k[z? + ?]'”) sin (20/3) dt = O(1) 


o 


for all z and all y < 0. 
On the other hand, 


_— ; - ) r,2 9° so 
v*(x, y) = e* [ ~*§ —_. K (ele? 4 271/2) 4 
ai Y) be € ot ox o\ | TT Zz ] ) Cc 
a 0 ke  s 2341/2 v ” =f rs) = 2 2941/2 
a Kolk[a + yy)” ) +. e —K(klt + 2°]”) dt 
ox Ja Ox 


fa} < — 
— K,(kr) + O(]) 
Ox 


Il 
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forx + y°—0,y < 0,as can be seen by the argument used to evaluate L, in the previous 


section. Hence 


e / e~'s*(t, 9) dt = Be? / a. 3 K(kié + y?}) dt + 0) 


Ko(kr) +e [ “e'Ka(kl® + y?}'”) at + OU) 


7@ 


K (kr) + O(1) 


forz + y — 0. Hence we obtain 


ole, v) = = BiKolkr) + O(1) 


for r — 0. Thus if B, = O we obtain a solution which is finite at the origin; if B, + 0 
the solution will have a logarithmic singularity at the origin. 

To complete the verification of the condition (2.2) one uses (5.11) to estimate the 
derivatives. 

In summary, for k < 1 there are two linearly independent solutions satisfying the 
conditions (2.1) and (2.2). One can be taken to be finite at the origin. The other then 
has a logarithmic singularity there. At infinity for y = 0 both of these solutions have a 
sinusoidal form. On the other hand for k > 1, there is only one linearly independent 
solution satisfying the conditions (2.1) and (2.2) and it has a logarithmic singularity at 
the origin. For k = 1 with y = 0 this solution approaches a constant for zr — + @; 
for k > 1 the solution dies out exponentially as x — . In every case the solutions die 


out exponentially as y — — ©. 
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QUASI-TRIDIAGONAL MATRICES AND TYPE-INSENSITIVE 
DIFFERENCE EQUATIONS* 


BY 
SAMUEL SCHECHTER 
Institute of Mathematical Sciences, New York University 


1. Introduction. In solving linear partial differential equations by finite difference 
methods a boundary value problem is reduced to solving a set of linear equations. 
In such instances the matrix involved usually takes a special form and consists mainly 
of zeros. Many of these matrices fall into the class to be considered here which may 
be called quasi-tridiagonal matrices. That is, we consider partitioned matrices of the 


h 





lorm 
[M, BE, O «+ =.+ - 0 | 
1D, M, E, 0 
Q=|0 D, M, FE, 0 + + O}=[D.,M.,B)f, (10 
0 . . - © De, Mey Bur 
0 . . ‘os 0 D, M,4 


where the D, , VM, , Z, are matrices with the same number of rows, EF, , M,.; , Da+2 
have the same number of columns, and the 1/7, are square. We propose to solve 


Qv=q (1.2) 


by direct methods. 

Various discretizations lead to matrices of this type. That the usual finite-difference 
approximation of certain boundary problems for the Poisson and the bi-harmonic 
equation yield matrices of the form (1.1) has been shown by O. Karlqvist [4], A. F. 
Cornock [2], L. H. Thomas [8] and others [6, 7]. Whereas the usual method of solution 
of the resulting linear equation is by iteration these authors propose direct methods 
for solving the above problems. 

The processes described here are, in part, extensions of those described by Karlqvist 
[4] and Cornock [2]. Furthermore it is shown that the finite difference equations obtained 
from symmetric positive systems, as defined by K. O. Friedrichs [5], also fall into 
the class of matrices of the form (1.1). These include, in addition to pure elliptic or 
hyperbolic equations, a certain class of boundary problems for equations of mixed 
type such as the Tricomi equation. Where iterative methods for such problems seem 
to present difficulties, even when Q is positive definite and symmetric, the direct methods 
for solving (1.2) are shown below to be feasible for this larger class of problems. 

A criterion is given for the process to apply which is similar to that found for the 
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tute of Mathematical Sciences, New York University, under Contract AT(30-1)-1480 with the U. S. 


Atomic Energy Commission. 
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LDU theorem [3]. It is also shown that when the MW, have the same order p, and the 
D,, are easily invertible, then the process may be reduced to multiplication of matrices 
and the inversion of one matrix of order p. This last fact was noted by Cornock [2] 
for the Poisson and bi-harmonic case. 

More generally, if forr = 1,2, --- ,k;k < q, and for integers q, such that 1 < g, < 
G2 <<*** <Q = 4, % = O, the matrices M, , for g,-1 < n < q, , all have the same 
order p, then it is shown that the process may be reduced to multiplications and the 
inversion of k matrices of orders p; , «++ , PD: . 

A code to solve (1.2) has been written by Max Goldstein for the I.B.M.-704 at 
New York University and this code has been successfully applied to a symmetric positive 
problem for the Tricomi equation. 

This code has also been applied to solving pure elliptic problems to compare the 
direct method with iterative methods. The direct method was used, for instance to 
solve the bi-harmonic boundary problem of a simply supported rectangular slab. A 
comparison of running time would indicate that the direct method is considerably 
faster than iterative methods for this type of problem. 

2. Direct methods. We seek a reduction of Q to the form 


Q = LU, (2.1) 

where Z and U are square matrices, partitioned in the same manner as Q, of the form 
L={C.,f.,05 , (2.2) 

U = (0, A, ,£, I, (2.3) 


t 


where J, is a unit matrix of the same order as MW, . We also partition the column vectors 


v and g in the same manner. 


Vv; Ji 
Vo | Jo | 
=| ge) ” 
. | 
Va L9e 


Comparing right and left hand sides of (2.1) we set, for 1 <n < q, 


A, = M, ? D, ”- Cohen ] M, - CLE .-1 + A, ° 


If the A, are non-singular, the A, and C, may be obtained recursively from 


A, = M, — D,A,-E,-: ; A, = M, (2.4) 
Co = Bs l<n<q. (2.5) 
To solve for v let Uv = y, then Ly = g; y and v may then be obtained recursively from 
Yn = On — CrYn-1 : tn = @, (2.6) 
where 
Y= ji ; and 


v, = Alyn — Endo Li" < @, (2.7) 





1960] QUASI-TRIDIAGONAL MATRICES 287 


where 
a ers 
Ue = Ag Yc: 

We will refer to this recursive method (2.4)-(2.7) as an LU-process. It should be 
noted that this process requires the inversion of A, , A, , --- , A, and their storage 
for use on the “backward sweep”’ (2.7). 

In certain cases most of these inversions may be avoided. For instance, if all the 
matrices M, are of the same order p and the D, are easily invertible, we may premul- 
tiply Q by the quasi-diagonal matrix 

I 0 | 
| Dp,’ | 
2 


| D3" 
| ° | 


| 
Lo D;'| 
We may thus assume that D, = I, 2 < n < gq. If we let A, = H,2,H, then the LU- 


process becomes 


H, = H,..M, — d...f.-; 2sn<Z 4, (2.8) 
where H, = 1,H, = M;; 
Hi-1(Yn — Gn) = —Ha-2Yn-1 5 2<n<q, (2.9) 
where y; = g, and Hv, = H,-,y,. 
If we let u, = H,-_,y, then from (2.9) 
Uy = HyisJn — Un-i » l<n<4q, (2.10) 
where u,; = g; . To obtain the v, we go back to the original equation (1.2) whence 
Un-1 = gn — Myv, — E,nvass ; li<n<q-1 (2.11) 
where v,_; = g, — M,v, and », is solved from 
Hw, =u, - (2.12) 


We have, in this case, only one matrix H, of order p to invert. 

This recursion (2.8)-(2.12) will be called an H-process. This method may be used, 
for instance, in solving problems for the Poisson or bi-harmonic equations over a rectangle. 
In the case of Poisson’s equation the qg inversions of the LU-process are reduced to 
solving only one set of p equations where p is the number of mesh points on a horizontal 
line. This fact was noted for these equations by Cornock [2] in specific examples. 

For regions which are made up of rectangles a similar result may be obtained. For 
instance, in solving Poisson’s equation for an L-shaped region made up of two rectangles 
R, , R, each having p,; , p2 points on a horizontal mesh line, respectively, only two 
matrices of order p, and p, need be inverted. 

To treat the general case let us assume that for 1 < q, <q@2< ++: <<Q@=9,% = 0 
the M,, have the same order p, for q,-; < n < q, . The matrix Q can then be partitioned 
again by combining those M, having the same order. That is, let 
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Q = (Di, Mi, Em, 





where 
[TO es) O Dass] 
D'., = 0 | l<r<k, 
LO 0 J 
Mrar=(D,,M,,#,)071, OSr<k 
0 0 
| O0O<r<k-1 
0 


LE..., O «e+ OJ 


, 


The LU-process can be performed for this new form for Q. Let g/ ,y/ ,v},l1<rsk 
be corresponding column vectors of dimension p, repartitioned as Q, where, for instance, 


we denote 


Var+1 
v1 = | ere a 
i 
Let yf = gf, ys = 92 — Chy{ where Dj = CiA{ . If we let 
yi = Alw! (2.13) 
w{ may be obtained by an H-process for (2.13). If we denote the H matrices that enter 
here by H, , --- , H,, then only H,, need be inverted and 
[ ges1 — De, a1We, | 
Y3 = | Ja.+2 - 
Lo? J 


Thus it is seen that only the last vector component w,, of w/ is needed and the reverse 
sweep given by (2.11) may be omitted. It will also appear later that only w,, and H;'H,,-: 
need be saved for future use. (It should be noted that where an H process is used the 
D matrices in M/ are assumed to be non-singular and that this process is feasible only 
if these D’s are easily, if not explicitly, invertible.) 

To obtain yj we must compute Aj , which would generally require the inversion 
of A/ . This however is not the case here. Since if we decompose A/ by the method of 
(2.1), (2.2), (2.3) into Af = L,U, then a straightforward computation shows that 


ge oe wo Bees 5 


DADE = DiU;'L EK! = |" 
_0 0 --- O 
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Thus Aj differs from M{ only in that W,,., is replaced by 

M,,+1 - D,, + HA ,-1 Ses ° 
That is Aj is also of the form (1.1) and has all of its M matrices of the same order p, . 
We are therefore reduced to the previous case. We may proceed in this manner until 
all the y’ , 1 < r < k have been computed and we note that v{ = wi . In this last com- 


putation, of course, (2.11) will have to be used to get all of vj . 
To obtain the other v! we note from (2.7) that, for r < k, 
vo) = wi — (A!) "E,, = wl — 21, 
where we let z! be the solution of 
Alz) = Ets, . (2.14) 
If the last vector component v,, of v! would be known the other components would 
be obtainable from the original equations as in (2.11). Since 
0 


0 

Evi. = | : 
0 | 
Rete cal 

if follows from (2.10) and (2.12) that if an H-process were carried out for (2.14) then 


-1 7 
4 = Hyg, -E eens: ? 


“ar 





where H;,'H,,-, is the matrix that entered in the computation of w,, . Thus if these 
are saved from before, the remaining v/ are obtained without any further matrix in- 
versions. 

Assuming therefore that the work to invert the D-matrices in the M{ is negligible 
we see that for a problem involving k sets of 1/-matrices of equal order only the k matrices 
H,, , H,, , +++ , H, need be inverted. Thus in solving Poisson’s or the bi-harmonic 
equation over a region made up of k adjoining rectangles the inversion of at most k 
relatively small matrices is required. A similar statement applies in higher dimensions. 

We note that in the particular case where Q is the discrete Laplacian, where M, = J, 
n= 1,2,--:,9,D,, E, are identities and J = [1, —4, 1]? , the inverse of H, can be 
given explicitly [4]. If we denote the allied Chebyshev polynomials by 


sinh (q + 1)z 
sinh x 


h(a) = , 2 cosh x = a, 


then H, = h,(J). The eigenvalues of J and H, are given by 


An = 2 cos ee) — 4, hom); m=1,2,--+ ,p, 


respectively, and the matrix of normalized eigenvectors is 


2 ms kmr y 
G = (2) ine PO 
a MO waa ey | 
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If we let 
L= [h.(A1) ee » ha,)] 
be a diagonal matrix then 
Ho’ = GL""G. 


For large problems the H, may be ill-conditioned. For the above problem the P- 
condition P(H,), or ratio of largest to smallest eigenvalue of H, is 


_msinh 6 (¢ + 1) 








P(H,) ~ + ; 
(p + 1) sinh 6-sinh (24 r) 

so that this may get quite large for large g = p. The inversion of H, may then present 
severe difficulties. It may however be feasible to break the problem into & groups as 
indicated above even though all the M, have the same order. 

Since (2.11) represents a marching process there is also the possibility of severe 
loss in accuracy in the value of v, for large g. The value of p does not appear to be an 
important factor in this loss for the discrete Laplacian problem. In the cases tried for 
gq = 5, a 704 code yielded results accurate to at least five significant digits for values 
of p= 5, 10, 20, 40. 

3. Criterion for decomposition. A sufficient condition for the validity of the de- 
composition is similar to that given for the LDU theorem [3]. Let 

M,E 
Oo, = M dQ. = ( or 
21 | 22 D.M, , ’ 
+ 1k 
Q, = [D,, ? M, » EL, ); eo Q, = Q. 

Tf Q: , Qe , +++ 5 Qa are non-singular then the decomposition (2.1)-(2.3) exists, and 
is uniquely given by the recursion (2.4), (2.5). In this case det Q = Bit det A, and if 
the M, all have the same order and D, , --+ , D, are non-singular then det Q = det H,. 

The proof follows easily by induction. From the Schur-Frobenius formula 


det G..; = det bs K | 
R Me, 


= det Q, -det (Mis, — RQ;'K) 


[ 0] 
| 
where R = (0::+0D,,:), K =| - |- From the inductive hypothesis 
| | 
ae 
LE, 
| 
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so that RQ7*K = D,,,A7’E, and det Q,., = det Q,-det A,,, . This implies that A,., 
is non-singular and also proves the formulas for det Q. 
A sufficient condition for Q, , --- , Q, to be non-singular is, clearly, that for all u 


yu'Qu > u"u (3.1) 


for some non-zero constant y. 
4. Applications. We consider, as an application of the above processes, the problem 
of solving a certain boundary problem for the Tricomi equation’ 


Tu = Yuz, — Uy — f(x,y) = 0 (4.1) 


by a difference approximation given by Friedrichs [5]. These difference approximations 
for symmetric positive systems have been further investigated by C. K. Chu [1]. 
The problem for (4.1) is posed for the parallelogram 


P:|y—2|<t, lz|<r;t, r>0O (4.2) 
such that 
Tu = 0, (x, y)eP (4.3) 
u,tu,=0 for |xr—y|=t, lei<¢r (4.4) 
u,=0 for x= -1, ly—ax|<t. (4.5) 


No condition is specified on x = r. 
Since the treatment given by Friedrichs calls for a rectangular region, P is trans- 
formed by — = x, n = y — 2 into the rectangle 


lé|<r, Ialse (4.6) 


and the equation (4.1) is written as the system 


" hs - 4 1 = 4 . (4.7) 
0 1 ; 63 0 


where y = £ + n, v = (v), 02)", 0) = Uz, Ve = U,. It is shown in [1] and [5] that after 
a premultiplication of (4.7) by a suitable two by two matrix of the form 


es 
1 p 


(4.7) and its boundary conditions can be brought into the required symmetric positive 
form. That is (4.1), (4.4), (4.5) can be written in the form 


d(a‘v, + a'v, + (a'r); + (a'v),) + xv = g, (4.8) 
where for a constant e 
Su m | eo (4.9) 
y p 





1A subscript of z, y, = or n indicates partial derivative. 
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is ee ¢ (4.10) 


pry Te p} 
of 
g= |", (4.11) 
f 
<n ia 9l' TO ‘|. (4.12) 


\€ €) 


As is shown in [1] and [5] the boundary conditions can be written in the form 


4 


Bv = pv, (4.13) 
where 

[ 
8 = —a’, B= -al=- py y | for €= —r, (4.14) 

2y 

| aes 

p 

B =a‘, u=ati =a* for t=r (4.15) 
B= —c’, w= —a. for » = —t where (4.16) 


— 1 (207 _ y(p +1) (l+ p)(p = 4 
é. 2 ee | 
[d+ p(e-y) p—2y+1) 


i] 


B=a’, p=a! for 7=t. (4.17) 


p-l 


To obtain the difference equations we divide the intervals (—r, r), (—#, ¢) by the 
2p — 1, 2g — 1 points, 








&; = (i — p) dé, 1l<4+sD— i], : € 9, (4.18) 
n=(G-9 An, l<js2q-1, 1<4q, (4.19) 
respectively, where p and q are odd integers and 
r t 


are the respective interval lengths. 

According to [5] an equation is written only for the pg odd numbered points (£; , 7;), 
i=2m—1,j=2n—11<m<p,1<n< q. If we write v;; = v(&; , 0;), af; = 
a®(t; , n;) and similarly for the other variables then the difference equations are given 
for an interior point 1 << m <p,1 <n <q, by 
vr eer , | 
— (0541, 542.5 — Oi-1,%s-2,5) + Ok (a ¢4:1Vig42 — OF,7-105,.7-2) 1 Kisii; = Giz, (4-20) 
where h = 2 Aé, k = 2 Ayn. 

For a boundary point at least one of the subscripts 7 + 1, j + 1 falls outside the 
range prescribed. In that case the a and v of the corresponding term are both evaluated 
at that boundary point and then replaced according to the rule (4.13). The 2h and/or 
the 2k in the difference involved are also replaced by A and/or k, respectively. 
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As an illustration we consider the equation for the point (£, , ;), not a corner: 
] , £ r_&) 1 n a] 
h (a9, Us, j = (a: )sitis) + a (at, 54101442 o 1, j-1V1, ;-2) + Ky ;01; ae 9ii ° (4.21) 


For the corner point (& , m), j, 2k, a], ;-:01,;-2 are replaced by 1, k and (a2)i.01, , 
respectively in (4.21). The other types of boundary points are treated in a similar manner. 
Since each (m, n) yields one equation we obtain pg equations in the pq unknowns 
»,; , or 2pq scalar equations for 2pg scalar unknowns. For each of the equations we 
note that v;; appears with at most four of its neighbors v,42,; , Ui-2,;  Ui,;-2 Ui.i+2 » If 
and denote by v” the vector with 2p scalar components arising from the 


we set v5 = 0;; 
;) on a horizontal line, the difference equations take the form 


points (&; , Non 





dnvm | + anvn-1 + div + Cvs + ene = Gn (4.22) 
ism Sp, l<n<q. 
This system can be written in the form (1.1), (1.2) where 
M, = [a., , bn, Calmear » 
D, = [0,d, , O}n-1 , 
E,, - [0, en . a ’ 
"| g’ | 
| y? 2 
v= | = 9 | ’ 
& | 
bs 
Lv | Lg’ 
and where forl1 <<m<p,l<n<q 
1 n l . 
d,, = ny Renee Om = 5p Mi-1.3 
n n 1 E 
by. = K; ’ Cm _ an 1,2 
Pe 1 ” rn 
m = 3), @: i+ Ym = iy 
2k ° & »® g g 
For the boundary points n = 1, g, 1 < m < p, 
di, _ 0, d,, = —* Qi ,2¢-2 
1 rn ¢ 
Cm = EM?» e, = 0 
andform=l1,plsn<q 
n n 1 £ 
= 0, a, = ~ fy S20- 23 
1 " 
q = — as ’ Cy —_ 0 
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Form=1,pl<n<4q, 
l ’ 1 
b, = k ZO)ii se On = hig HF owes 
and for corner points 
Pe ‘ee 
b Ki3 — Tt (a)as — = (a_)3; | m=] n=1 
h NO el P , 
rae 7 
b Kg bP aii — 7 mM =D, n= 1 
a es OF 
. _— h (a-)i3 T ] \Q-)ij 5 m=1, eS hy 
. cs . §.. 
Og Ky or = Gis oT os ei 5 m=D, n= q; 
i h ] k I 
where in all cases? = 2m — 1,7 = 2n — 1. 
Thus this boundary problem for the Tricomi equation yields a matrix of the form 


(1.1). It has been shown by Chu [1] that, if «, r, ¢ are properly chosen, e.g., « = 1/2, 
r = 1/2,¢ = 1/5, this matrix satisfies the inequality (3.1). The LU-process is therefore 
applicable. Since the M, all have the same order we may in fact use the H-process 
providing the D, are non-singular and easily invertible. The D, are quasi-diagonal, 
and a simple computation shows that for the above choice of e¢, r, ¢, a’ is non-singular. 


The above mentioned code was used to solve the problem described above for the 
choice of 
f(a, y) = 6yx — 4y° + y —1 — 2x 
and was run for various values of p and g. The solution to the analytic problem (4.3)- 
(4.5) can be given explicitly by 
u(z, y) = (y — 2° (4 +2) — 2/25 (4.23) 
and for (4.8)-(4.17) it is given by 
. ia 2 a ( 1 "eo E 
vi n nl + 2€) 1/25, (4.24) 


“aa3 


v2 = n(l + 2é). 
For the case of p = 15, g = 15 the code yielded an answer accurate generally to 
three significant digits. The values given at J(0, 0), JJ(1/2, 0) by the code are illustrated 


by Tables I, II respectively: 











pxq __ 3 xX 3 5x5 7 X13 11 X7 15 X 15 es 
Y, —6.22 —4.07 —4,.02 —4,.27 —4.04 
I x 10° 
ve | ———0.37_ 0.08 0.01 —0.04 —0.01 —_ 
v; —7.13 —4.96 —4.34 —3.99 —3.99 
Il x 10°? 
V> —1.01 0.54 0.14 —0.21 —0.05 


The exact solution given by (4.24) for 7 = Oisv, = —.04, v. = 0. 
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Other problems for the Tricomi equation were run to check the influence of the 
positivity and boundary conditions. In one case the two by two premultiplier matrix, 
needed to guarantee “positivity,’”’ was omitted. Again an approximate solution was 
obtained with a slight loss in accuracy. For instance the value at J given by the code 
was ?; = — .0454, vo = — .0019, for p = i, qg= a 

A problem for the homogeneous Tricomi equation with inhomogeneous boundary 
conditions (u — 8)v = f, for a given f, also yielded results similar to those given above. 

Symmetric positive systems for dimensions higher than two may be treated in a 
similar manner. 1 

As is pointed out in [2] and ry, problems for the bi-harmonic equation can also be 
put into the form (1.1). The LU-method was carried out for a simply supported rec- 
tangular plate and was found to give results accurate to about three significant digits 
for a 30 by 30 mesh. The running time for this problem was about an hour on the IBM-704. 

The above methods may also be combined with iterative or group relaxation methods 
where each individual group relaxation is done by a direct method. This has already 
been proposed for a multigroup diffusion problem by Nohel and Timlake [6]. 

It may also be noted that, in problems where higher order difference schemes are 
available, the direct methods will require relatively little additional operations and 
thus one may require fewer mesh points in a given problem. 
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—NOTES— 


A NOTE ON RIGID BODY DISPLACEMENTS IN THE THEORY 
OF THIN ELASTIC SHELLS* 


By P. M. NAGHDI (University of California, Berkeley) 


1. Introduction. The stress-strain relations appropriate to the linear theory of thin 
elastic isotropic shells have long been the subject of numerous investigations. Recently, 
a suitable system of stress-strain relations, which include the effects of both transverse 
shear deformation and transverse normal stress, has been derived for axisymmetric 
shells of revolution by E. Reissner [1] and, in general form by the present author [2], 
where the deformation is referred to the lines of curvature of the middle surface. Also 
contained in [2] are explicit expressions for the stress resultants {N, , N., Ni2, Na}, 
{V, , V2}, and the stress couples {M, , M. , My. , M2,} in terms of quantities which 
describe the state of strain, when only the effect of transverse shear deformation is 
retained (see Eqs. (3.6) of [2]). This latter system of stress-strain relations, in the notation 


of Ref. [2], is of the form 
T y 0 0 h? 
N 3 >= an Y2) + — (et -; aN R, = i) |, 
0 0 ] (2 
N., = = onl +7) += “(k- R, 
l—~yp ] ] 
M1. = sae D (3 + 52) + (7. - 1)], 
2 R, R, (1.1) 


1 
My, = -—— 5 Dl va, + 8) + (4 -% he], ete., 


Cn” 


which satisfy the identity 


ee 
Ni + 2 = WN; 


In (1.1), R, and R, are the principal radii of curvature of the middle surface; h is the 
shell thickness; D = Eh*/[{12 (1 — v*)]; G is the shear modulus; EF and » are respectively 
Young’s modulus and Poisson’s ratio; and the quantities y? 6, , etc., which refer to the 








components of strain, will be defined presently. 

The stress-strain relations derived in [2] have been criticized very recently (by 
Goldenveiser [3], who has asserted that for rigid body displacements alone not all quanti- 
ties on the right-hand side of (1.1) vanish, which implies that the shell space will not 
remain stress-free. It is the purpose of the present note to examine the implications of 


rigid body displacements of the shell space, and this may be regarded as an addition to 





*Received March 19, 1959; revised manuscript received January 6, 1960. The results presented in 
this paper were obtained in the course of research sponsored by the Office of Naval Research under 


Contract Nonr-222(69), Project NR-064-436, with the University of California (Berkeley). 
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the results in [2]. Although attention is confined to shells whose stress-strain relations 
are defined by (1.1), the extension to the case in which the effects of both transverse 
shear deformation and normal stress are included will be apparent and leads to the 
same conclusion. It is shown that in the absence of surface loads and body forces, con- 
trary to the assertion made in [3], for rigid body displacements all of the stress resultants 
and the stress couples vanish; and conversely, the vanishing of these resultants corre- 
sponds in fact to the most general form of the rigid body displacements of the shell space. 

2. Preliminary background. In order to make the following analysis reasonably 
self-contained, we recall that the position vector of a generic point in the shell space 
is defined by 


R¢é, , & , $) = r(é, , &) + gn (2.1) 


where r is the position vector of the corresponding point (with coordinates £, and &) 
on the middle surface, and — h/2 S ¢ S h/2 is measured along the unit normal n to 
the middle surface. Further, we denote by t, and t, (t, , t, and n form a right-handed 
system), the unit tangent vectors to the é,- and £-parametric curves which will be 
referred to the lines of curvature on the middle surface (¢ = 0). 

The displacement vector, as assumed in [2] and which leads to the results (1.1), is 
given by 

VE; ’ £> ’ 9) — u(é, » €2) + cB, , £2), (2.2a) 

where 


u = ut, + Ut» + un, 8 = Bt, + Bot, . (2.2b) 


With the notation 0/d&, ( ) = (_ ), , the quantities y? , 6, , ete. in (1.1), (or in Eqs. 
(3.6) of [2]), may be conveniently written as 





0 u 1 J 1 0 u 2 r 2 
5 = ae 3 = —-— © 
QQ) Q2 
5. Bats 0 West, 
a= m= 
QA QoQ, 
ayn: = Br, + Un, QsxYor = BT.» + U.-n, (2.3) 
6.4°T., Y.o'D 2 
is ¥ ’ a=” "3 , 
Qa, Qo 
¢ 8 iT 2 8 2°T 1 
Qs e'. eS 4 
a)Qe AQ; 
where (a? , a3) are the non-vanishing components of the middle surface metric. 


3. Rigid body displacements of the shell space. Let, in the absence of surface loads 
and body forces, the displacements be specified by 


U=u+0XR, (3.1a) 
which, by (2.1), may also be written as 
U=(u+ Xn + fo Xo), (3.1b) 


where uy and the rotation » are constant vectors. The displacement vector defined by 
(3.1) represents the most general form of rigid body displacements of the shell space, 
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to 
oo 


and includes as special cases purely rigid body translation and purely rigid body rotation. 
Yeturning now to our main objective, we substitute (3.1) into the first six of (2.4) 
to obtain 


’ 1 
é— = “ef 3°(@ x r ;) = 0, 
ay 


"1 = —12‘(o XPF,;) = on, 
Qa;:ao 
saptie 3.2a) 
] if , ! 
YY =— fio Xn+n-o Xr]; 
ay 
es 
=—[|r,-o Xn—r1T,-o Xn] = 0, 
Qa, 
and similarly, 
e=0, y,=--en, yx = 0. 3.2b) 


Likewise, substitution of (3.1) into the remaining expressions of (2.4), with the aid of 
N,, = fF, R, and ni. = T,o i 9 yields 


l l 
=—-30 X,T, = zen, Xr, = 0, 
Qh) ay 
, ] l ee re 
i =- OA Bits = 6 - s 3.3a) 
Q1Q> Q)Q2 R, 
a on 
im ¢ 
and similarly, 
a o-n Q « 
Ko = U, = -—=-: 3.3b) 
R, 


Next, substituting (3.2) and (3.3) into (1.1) and observing that 


WwW 
— 


ee (J hye n+a+(2 ae 0 
(6, + 6) + (> —- pm = (64 + &) + 1S — oN: = 9, 
:, R, | pe dy R, R, U 
it follows that in the absence of surface loads and body forces, all of the resultants 
in (1.1) vanish. Furthermore, it may be easily verified that when in the absence of 
surface loads and body forces all resultants are set equal to zero, the solution of dif- 
ferential equations (1.1) and (2.3), with an appeal to the uniqueness theorem, is indeed 
given by (3.1). 
Note added in proof: Professor Goldenweiser has recently informed the author that 


he is in agreement with the content of this note. 
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A NOTE CONCERNING HYPERBOLIC EQUATIONS 
WITH CONSTANT COEFFICIENTS* 


By JAMES F. HEYDA (Aircraft Nuclear Propulsion Department, General Electric Co.) 


The canonical form of the linear, homogeneous hyperbolic partial differential equa- 
tion in two independent variables is 


Uzy + au, + bu, + cu = 0, (1) 


where a, b, c are functions of x, y. The constant coefficients case, frequently met in appli- 
tions, can be reduced to the simpler form 


v,, — vw = 0, A=ab-—c (2) 
with \ constant, through the substitution 
v(x, y) = u exp (bx + ay). (3) 


The solution of the Cauchy problem for an equivalent form of (2) with v and dv/dan 
prescribed along C is worked out explicitly by Copson [1] using the methods of M. Riesz. 











&,, 
(*., 4%) Y 








Fia. 1 


The form of Copson’s solution is the same as that which one would obtain by using 
Riemann’s method for solving the same problem. This enables us to identify the Bessel 
function J,(2[A(z — a)(y — 8)]'”’) as the Riemann characteristic function for Eq. (2) 
By a theorem of Le Roux [2], it then follows that the definite integrals 


*Received September 19, 1958; revised manuscript received August 28, 1959. 
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[. fla)To(2[My — yo)(x — a)]'*) da, 
Zo (4) 


[ s@tene - xy — at”) a8 


will also be solutions of (2) for arbitrary functions f(a) and g(@). This enables us to write 
down directly a solution of (2) for prescribed values of v(x, y) along the characteristic 
lines x = x) and y = yp . For we have from (4) by superposition 


az 


v(x, y) = | fla)I,(2[My — yo)(x — a)]'’*) da 


+ | g(B)I(2[Ma — ao)(y — B)]'”) dB + v(ao , yo)Lo(2[A(e — roy — y)]'” 


where 
du(Xo » y) 


f(z) = See gy) = —— 6 
LM, ax Jy ay ’ ) 





the prescribed functions v(x, yo), v(%o, y) being assumed suitably well-behaved to carry 
out the operations indicated. 

As an instance of (5) we remark that it furnishes directly a closed form solution of 
the boundary value problem treated by Mason [3] in connection with heat transfer in 
cross-flow and explains the mysterious appearance of J,(2[abzy]*) in his solution, arrived 
at by means of Laplace transform methods. 
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Corrections to the paper 


THERMAL INSTABILITY OF VISCOUS FLUIDS 
Quarterly of Applied Mathematics, XVII, 25-42 (1959) 
By CHIA-SHUN YIH (University of Michigan) 


Equation (27) should read 
[o — Pr (D* — b’ — a’)] (D* — a’)f = —R Pr DO 


The second sentence following Eq. (27) should then be changed to: “Since Eqs. (26) 
and (27) would be identical to Eqs. (20) and (21) if ¢ + 6° is replaced by a, any variation 
of the flow with y will invariably make the flow more stable, if uw, remains strictly zero, 
and the boundary conditions on f itself are not relaxed.” 

With the periodic variation in the y-direction provided, the boundary, conditions 
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on f itself can be relaxed, and the fluid may be more unstable. This instability in fact 
corresponds to the root zero of the characteristic equation 


tan R'* = tanh R™* 


given in the paper. The author is indebted to Mr. R. A. Wooding of Cavendish Laboratory 
of Cambridge for pointing out this most unstable mode. 


ON THE SOMMERFELD HALF-PLANE PROBLEM* 
By BERNARD A. LIPPMANN (Lawrence Radiation Laboratory, University of California, 


Livermore, California) 


Abstract. A simple derivation of the Sommerfeld solution to the problem of the 
diffraction of a plane, scalar wave by a half-plane is given. The discusson is of interest 
mainly because of the simplicity of the argument; however, it is felt that the ansatz 
that forms the core of the argument is probably more generally applicable. 

I. The method of deriving the Sommerfeld solution to the problem of the diffraction 
of a plane, scalar wave by a half-plane, presented here, became apparent during a 
study of the utility of conformal mapping in diffraction problems. The discussion is 
presented mainly because of the simplicity of the argument employed to deduce Sommer- 
feld’s well-known result, but, beyond this, there are two other features of interest: 
first, one feels that the ansatz that constitutes the core of the argument is probably 
applicable more generally; and second, although the explicit conformal mapping used 
is trivial, the manner in which the mapping enters into the formulation of the boundary 
conditions in the ansatz may be suggestive in clarifying the relationship between con- 
formal mapping and diffraction theory. 


II. We seek a function, u, that is a solution of the wave equation’ 


(2 _- = + ae = 0, (1) 


Ox” ay" 
corresponding to the incident wave’ 
exp [—7kr cos (6 — 6’)], (2) 
and satisfying certain other conditions to be stated presently. These conditions, the 


wave equation, and the form of the incident wave, show that u is a function of 6 — 6’ 
only; we may therefore allow 6’ to approach zero and use 
Uy = exp (—7kr cos 6) (3) 
as the incident wave, providing that, in the final result for u, we replace 6 by 6 — 0 
everywhere. 
It is instructive to consider the problem simultaneously as it appears in Fig. 1, 
which shows the actual half-plane, coincident with the positive x-axis, and Fig. 2, which 
*Received Oct. 26, 1959. Work done under the auspices of the U. S. Atomic Energy Commission. 
!We shall use the coordinates z, y; r, 0; and z = y + iy, z* = x — iy, as convenience dictates. 
*In contrast to the usual convention, note that this incident wave is traveling towards the z-axis, 


from above, at an angle 6’. 
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shows the half-plane mapped on the entire real axis. (Figure 1 is related to Fig. 2 by 


z= 2, .) 

We shall assume that 6’ has approached zero from positive values of 6’, so that the 
incident wave in Fig. 1 is regarded as coming on the half-plane from above, right, at 
grazing incidence, and moving in the negative z-direction; in Fig. 2, the same wave, 
suitably transformed, is incident on the entire real axis from above. 

Suppose the boundary condition on the half-plane is u« = 0; that is, u vanishes all 
along the real axis in the z, plane. This condition may be satisfied as follows. Assume 
we find a solution of the wave equation that has wu, as its incident wave and that satisfies 
any boundary condition on the 2,-axis. From this solution we subtract the function 
obtained by reflecting the incident wave in the z,-axis; the result is a solution of the 
wave equation that has the correct source above the x,-axis and reduces to zero on 
the 2,-axis. 

Thus, the only condition we need impose on u, other than that it is a solution of 
the wave equation, is that the wave wp is incident on the half-plane from above; the 
boundary condition on the half-plane may be chosen at our convenience. 

One simple choice regards the half-plane as “black”. That is, the half-plane absorbs 
the incident wave completely. No scattered wave is required, except at the edge of the 
half-plane; here, a scattered wave is needed in order to maintain the continuity of the 
wave function u in going around the half-plane. 

This boundary condition means that we must find a (continuous) function u such 
that: 

u satisfies the wave equation. (A-1) 
u—>U, as xz,— © fory, > 0. (A-2) 
u—Oasz,— —-@ fory, > 0. (A-3) 
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We now proceed to satisfy these conditions. Our method is based on the fundamental 


ansatz: 


u = ul (z,), (4) 
where we require that 
Ko) =1, (5) 
and 
I(—) = 0, (6) 
a pair of conditions that can be met by putting 
I(x.) = | f(r) dr, (7) 


where a is an arbitrary constant. 
(A-2) and (A-3) are now satisfied. Moreover, f(r) is still arbitrary, leaving the pos- 
sibility that it can be chosen to satisfy (A-1). 





Since 
ee ae a 
ax” " ay’ | 2, |? de, det’ 
A-1) demands that 
1 9° ) ik a 
a 2 Ww aia ae a” — (+ 2* = 
>, |? de, oz +k bexp| 9 1 + 2 Ke ( +2 )) 0. (8) 
Irom this, we deduce that ie == @ 1 | 
aT 
f'(ax,) _ Atk mae 
flaxz,) a’ 
or, 
f(r) = @ exp (i2k7’/a’). (9) 


Putting this in (7), adjusting a to satisfy the normalization requirement (5), and using 
p”” cos 6/2, we obtain: 


(2kp)*/%e086/2 


nap'/2eo0s6/2 q 1/2 
I(x,) = | (24) exp (2ikr’/a’) = (wi)"'”” [ exp (i7’) dt. (10) 


rat Jaw 
The complete solution is now given by making the replacement @ — 6 — 6’ in (3), (4), 
and (10): 
exp (ir’) dr. (11) 


(2kp)*/2cos(0-—0')/2 
u(@; 6’) =(ni)~’’~ exp [—ikr cos (@ — 6’)] [ 
The solution corresponds to a “black” half-plane, but, following the prescription 
given earlier, it can easily be transcribed into the solution for a perfectly conducting 
half-plane (wu = O on the half-plane): In the z,-plane, reflecting the source in the real 
axis changes 6/ to — 6 ; this transformation replaces 6’ by — 6’ in the z-plane. Therefore, 

the required function is 
u(@; 6’) — u(@; — 6’), (12) 


since this takes on the value zero for 6 = 0; 27. 
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Finite-Difference Methods 
for Partial Differential Equations 


By Gerorce E. Forsytue, Stanford University, and WOLFGANG R. 
Wasow, University of Wisconsin. Of the many numerical methods for 
solving partial differential equations, only one stands out as being 
universally applicable to both linear and nonlinear problems—the 
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Markov process as a system model. 1960. 
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Classical Mechanics 
Second Edition 
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